Techniques for controlling charging and/or discharging of batteries using a tanks-in-series model

ABSTRACT

In some embodiments, a battery management system is provided. The battery management system comprises a connector for electrically coupling a battery to the battery management system, at least one sensor configured to detect a battery state, a programmable chip configured to control at least one of charging and discharging of the battery, and a controller device. The controller device is configured to receive at least one battery state from the at least one sensor; provide the at least one battery state as input to a tanks-in-series model that represents the battery; and provide at least one output of the tanks-in-series model to the programmable chip for controlling at least one of charging and discharging of the battery.

CROSS-REFERENCE TO RELATED APPLICATION

This application claims the benefit of Provisional Application No.62/905,942, filed Sep. 25, 2019, the entire disclosure of which ishereby incorporated by reference for all purposes.

STATEMENT OF GOVERNMENT LICENSE RIGHTS

This invention was made with government support under Grant No.DE-AR0000275, awarded by the U.S. Department of Energy. The governmenthas certain rights in the invention.

BACKGROUND

Lithium-ion batteries are achieving increasing ubiquity with theemphasis on renewable energy and electric transportation to addressgreenhouse gas pollution. Material and manufacturing innovations,coupled with market pressures and economies of scale, have led todeclines in battery pack costs of nearly 85% in the period 2010-18. Animportant factor towards enabling further reductions in cost metrics isthe efficient operation of battery systems. This entails maximizingutilization, reducing overdesign, enabling fast charging, and mitigatingdegradation phenomena. Thermal runaway and the resulting likelihood offires and explosions is another critical consideration, and it isessential to ensure battery operation in the safe temperature range. TheBattery Management System (BMS) is the combination of hardware andsoftware components that performs the requisite functions to ensure safeand efficient operation. An accurate and robust mathematical model isessential for BMS, which estimates battery state variables such as Stateof Charge (SoC), State of Health (SoH), and temperature.

The incorporation of sophisticated electrochemical models has thepotential to enable more powerful and intelligent BMS. In addition toimproved prediction of battery states, variables internal to batteryelectrochemistry allow the formulation of more complete optimization andcontrol problems than would be possible by simple circuit models.Achieving this integration requires the reduction of the computationaldemands of complex models, such as the macro-homogeneous p2D models ofNewman and co-workers. In addition to model reformulation techniquesthat exploit the underlying mathematical structure of the equations toachieve fast simulation, the tradeoff between accuracy and computationalefficiency has spurred active investigation into simplifications of thep2D model subject to limiting assumptions.

Since its introduction by Atlung et al., the Single Particle Model (SPM)has been extensively used for efficient simulation, estimation, optimalcharging and life-cycle modeling. The SPM visualizes the electrodes in alithium-ion battery as two representative spherical particles andconsiders the electrode reactions and transport through the activeparticle. SPM is particularly attractive since the largest reduction inthe number of equations is achieved via simplifications of the solidphase, given that the corresponding equation is solved at each point inthe electrode computational grid. Importantly, however, it neglectspotential and concentration variations in the electrolyte phase, whichlimits its predictive ability to low-current scenarios and relativelythin electrodes, where liquid phase polarizations can be neglected.

In addition, SPM assumes a uniform reaction distribution throughout eachelectrode, which is only attained at moderate currents in which kineticresistances dominate ohmic effects, and in electrodes which possess amonotonic dependence of equilibrium potential with degree of lithiation.SPM has been recovered in the limits of fast diffusion dynamics withrespect to characteristic discharge time.

The p2D model also returns SPM in the limit of large changes incharacteristic overpotential upon Li intercalation, which results inuniform reaction distributions. These limitations have led to efforts toexpand the applicability of SPM by introducing simplified descriptionsof electrolyte dynamics and non-uniform reaction distributions.

In past work by the inventors of the present disclosure, modelreformulation techniques have been used. Millisecond computation timeswere achieved using coordinate transformation, spectral methods andreformulation enabled by analytical solution, and successfullydemonstrated in applications for parameter estimation, optimal control,and BMS. Despite detailed analysis and comparisons of efficientsimulations of SPM-like models in past work, the two and three-parametermodels have been arguably the most widely used for control, optimizationand BMS by the community at large. This motivates the development ofefficient models for the electrolyte phase in lithium-ion batteries.

A common approach for the inclusion of non-uniform reactiondistributions is to directly assume a polynomial profile for pore-wallflux. Alternatively, polynomial profiles for both electrolyteconcentration and potential may be assumed and used to determine thespatial variation of the pore-wall flux. In addition, some workers havederived closed-form analytical solutions for the pore-wall fluxdistribution. However, these solutions are only valid under certainassumptions, such as linear kinetics, or by neglecting the diffusionalcontribution to the electrolyte current, and also neglect theconcentration-dependence of one or more electrolyte transportproperties. Still other approaches assume a polynomial or exponentialdependence for the equilibrium electrode potential in space. Theseassumed profiles are often combined with polynomial spatial dependencefor other electrochemical variables, converting the original PartialDifferential Equations (PDEs) of the p2D model into a system ofDifferential Algebraic Equations (DAEs). The simplifications useassumptions to achieve closure, while also assuming constant electrolytetransport properties for tractability. For both types of approaches, theuse of higher-order polynomial profiles for increased accuracy isexpected to increase the DAE system size, with the attendant penalty forcomputational efficiency.

Extensions of SPM to electrolyte dynamics are often termed ‘SPMe’, wheree denotes electrolyte. A common approach begins with the assumption of auniform reaction in the electrodes. This results in simplified PDEs forelectrolyte concentration. For constant electrolyte diffusivity, SPMeresults in linear PDEs, which are computationally simple to solve.Indeed, it is possible to derive closed-form analytical solutions forthe constant diffusivity case. Alternatively, instead of constantreaction distributions, polynomial spatial dependence for concentrationcan also be assumed, resulting in a system of DAEs. With the knowledgeof electrolyte concentrations, the electrolyte current equation isusually integrated, or volume-averaged quantities are used to obtain theliquid-phase potential drop. In some cases, polynomial profiles areassumed for the liquid-phase potential as well. This information is thenused in conjunction with the electrode-kinetics equation to estimate thecell voltage. The uniform-reaction SPMe is rigorously derived from thescaled p2D model in the limit of fast electrolyte diffusion dynamics.This perturbation approach yields simplified algebraic expressions forthe various polarization contributions to the cell voltage. Theseexpressions can be evaluated at negligible computational cost. Inaddition, the cell voltage equation is expressed in terms ofelectrode-averaged quantities to obviate the need for assumptions on thepore-wall flux at the terminals.

The PDE system in SPMe becomes non-linear in the case of variablediffusion coefficient, which might compromise computational efficiency.When solving numerically, the non-linearities are likely to entail finerspatial discretization in both the electrode and the solid particle toensure numerical convergence, which will increase the DAE system size.The stiffness of the resulting system is also likely to increase for thenon-linear case, especially at high current densities. The effect ofincreased computational cost may become more significant in real-timeenvironments, both due to stringent accuracy requirements and hardwarelimitations. In addition, even in situations where the non-linearequation is solved numerically, simplifications of the electrolytecurrent equation, such as constant electrolyte conductivity, arenecessary to allow tractable integration for calculating electrolyteohmic drop. This can be a requirement even if polynomial profiles forelectrolyte concentration are assumed for the non-linear equation. Theperturbation approach results in expressions in which the transportproperties are evaluated at a given representative concentration. Thismay result in errors in estimating concentration overpotentials at highdischarge rates, where substantial spatial variations are expected toarise. Additionally, it is not clear whether the algebraic expressionsfor electrolyte ohmic drop are valid for variable concentrationproperties and may have to be rederived. Using a constant value may leadto errors in estimating ohmic drops and concentration overpotentials athigh current-densities.

What is desired are techniques that can provide electrochemical batterymodels that can both be computed efficiently and can provide higheraccuracy and performance, in order to effectively and accurately controlcharging and/or discharging of batteries.

SUMMARY

This summary is provided to introduce a selection of concepts in asimplified form that are further described below in the DetailedDescription. This summary is not intended to identify key features ofthe claimed subject matter, nor is it intended to be used as an aid indetermining the scope of the claimed subject matter.

In some embodiments, a battery management system is provided. Thebattery management system comprises a connector for electricallycoupling a battery to the battery management system, at least one sensorconfigured to detect a battery state, a programmable chip configured tocontrol at least one of charging and discharging of the battery, and acontroller device. The controller device is configured to receive atleast one battery state from the at least one sensor; provide the atleast one battery state as input to a tanks-in-series model thatrepresents the battery; and provide at least one output of thetanks-in-series model to the programmable chip for controlling at leastone of charging and discharging of the battery.

In some embodiments, a method of managing at least one of charging anddischarging of a battery is provided. A controller device receives atleast one battery state of the battery from at least one sensor. Thecontroller device provides the at least one battery state as input to atanks-in-series model that represents the battery. The controller deviceuses at least one output of the tanks-in-series model to control the atleast one of charging and discharging of the battery.

In some embodiments, a computer-implemented method of deriving a modelfor controlling at least one of charging and discharging a battery isprovided. A computing device integrates pseudo 2-dimensional (P2D) modelrepresentations over volumes of a cathode domain of the battery, ananode domain of the battery, and a separator domain of the battery. Thecomputing device approximates interfacial fluxes between the domains,wherein a corresponding gradient is expressed as an average value minusan interfacial value of a dependent variable over a length scale tocreate a tanks-in-series representation of the battery. The computingdevice provides the tanks-in-series representation of the battery as themodel for controlling at least one of charging and discharging of thebattery.

DESCRIPTION OF THE DRAWINGS

The foregoing aspects and many of the attendant advantages of thisinvention will become more readily appreciated as the same become betterunderstood by reference to the following detailed description, whentaken in conjunction with the accompanying drawings, wherein:

FIG. 1 is a partially schematic view of a non-limiting exampleembodiment of a battery management system according to various aspectsof the present disclosure.

FIG. 2 is a block diagram of a non-limiting example embodiment of acontroller for a battery management system according to various aspectsof the present disclosure.

FIG. 3A and FIG. 3B are block diagrams that illustrate non-limitingexample embodiments of the present disclosure as they may be deployed indevices.

FIG. 4 illustrates the computational schematic of the pseudo2-dimensional (p2D) model of Newman and co-workers, a continuumelectrochemical model that has found substantial application forsimulation of Li-ion battery performance.

FIG. 5 illustrates a visualization of a non-limiting example embodimentof the Tank Model approach depicting the mass and charge flows in to andout of each ‘tank’ according to various aspects of the presentdisclosure.

FIG. 6 is a flowchart that illustrates a non-limiting example embodimentof a method of deriving a tanks-in-series model that represents abattery according to various aspects of the present disclosure.

FIG. 7 is a flowchart that illustrates a non-limiting example embodimentof a method of controlling at least one of charging and discharging of abattery according to various aspects of the present disclosure.

DETAILED DESCRIPTION

In some embodiments of the present disclosure, a Tanks-in-Seriesmodeling approach is used to reduce the p2D model of a battery. Theelectrolyte conservation equations are arranged in volume-averaged form,without resorting to any direct assumptions on the spatial dependence,unlike in polynomial profile methods. The proposed Tanks-in-Seriesapproach also does not assume uniform reaction rates to predictconcentration profiles, since we deal in average quantities. The keyapproximations are (a) in the interfacial fluxes and (b) using theelectrode-averaged pore-wall flux to estimate the electrode-averagedoverpotential. For problem closure, we make reasonable approximationsfor the flux variables at the interface. Mass and charge conservationare imposed at the domain interfaces in order to determine the unknowninterfacial values. The cell-voltage is then calculated from the knownelectrode-averaged pore-wall flux. This formulation allows for theinclusion of concentration-dependent transport properties, sinceterminal-to-terminal integration is not required, and the properties areonly evaluated at the domain interfaces. In addition, it reduces thefull p2D model into a fixed-size system of <20 DAEs and no PDEs need besolved.

An efficient, conservative model for lithium-ion batteries is presented,which uses a Tanks-in-Series approach to generate approximate equationsfor electrolyte dynamics. Despite the loss of information due tovolume-averaging, the additional equations result in a nearly four-foldreduction in error compared to typical SPM at high (5C) discharge rates.While the Tank Model achieves excellent accuracy with respect to thefull p2D model for the cell considered herein (<0.4% error), even betterperformance is expected for electrodes with inherently more uniformreaction distributions. The model formulation provides for theconvenient refinement of flux approximations and estimation to increaseaccuracy for more aggressive parameter combinations. The model retainsthe computational simplicity of SPM-like models, with the ˜millisecondcomputation time making it a candidate for a replacement of SPM in thesimulation of large series-parallel configurations of cells. The modelmay also be used for long cycle simulation and parameter estimationtowards modeling capacity degradation, and the evaluation of models forthe same. The generalized methodology described herein can be applied toelectrochemical models for more complex systems, including but notlimited to conversion chemistries such as Li-sulfur and lead-acidbatteries.

FIG. 1 is a partially schematic view of a non-limiting exampleembodiment of a battery management system 102 according to variousaspects of the present disclosure. The battery management system 102 isa non-limiting example of a system that may use the Tanks-in-Seriestechniques described herein to control charging and/or discharging of abattery.

In some embodiments, the battery management system 102 may include apower supply 104 and a DC supply 106 that provides power to anelectronics board 108. The electronics board 108 may include severalcomponents including a programmable chip 112 (e.g., an EPROM). A batterycharging housing 116 may include at least one rechargeable battery 114.In some embodiments, the battery charging housing 116 may be replaced byconnectors that are electrically connected to the electronics board 108with conductive wires.

The operation of the DC supply 106 and/or the electronics board 108 maybe controlled by a controller 110. For example, the controller 110 mayload data onto the programmable chip 112 that, in turn, controls the DCsupply 106. The data loaded by the controller 110 onto the programmablechip 112 may include one or more parameters for a model that describes,for example, charging/discharging, heating, cycling, etc., for therechargeable batteries. and/or may include the model itself. Forexample, parameters and/or models obtained by one or more of thetechniques described below may be used to control charging/dischargingcurrent, charging/discharging voltage, temperature of the battery, peakefficiency of the battery, optimal number of charging/dischargingcycles, etc., of the battery 114 through the controller 110 and/orprogrammable chip 112. In some embodiments, the controller 110 maycontrol the DC supply 106 directly. In some embodiments, the controller110 may be part of (e.g., may be carried by) the electronics board 108.

A single battery 114 is illustrated and described herein for the sake ofsimplicity. However, in some embodiments, the battery management system102 may support using more than one battery 114. If more than onebattery 114 is used, then there may be a separate controller 110 orprogrammable chip 112 associated with each battery 114, or a singlecontroller 110 or programmable chip 112 may manage all of the batteries114. One of ordinary skill in the art will recognize that the batterymanagement system 102 illustrated in FIG. 1 is an example only, and willrecognize that the techniques described herein may be used in othertypes of devices that are used to control charging and/or discharging ofa battery without departing from the scope of the present disclosure.

FIG. 2 is a block diagram of a non-limiting example embodiment of acontroller 202 (e.g., a computing device) for a battery managementsystem according to various aspects of the present disclosure. Thecontroller 202 may be suitable for use as a controller 110 asillustrated in FIG. 1 . The controller 202 includes one or more inputdevices 204 that provide input to a processor (CPU 210). Input devices204 can include, for example, a mouse, a keyboard, a touchscreen, aninfrared sensor, a touchpad, wearable input devices, a camera orimage-based input device, microphone, or other input devices. The CPU210 may be a single processing unit or multiple processing units in adevice or distributed across multiple devices. The CPU 210 may becoupled to other hardware devices, for example, with the use of a BUS,such as a PCI BUS or SCSI BUS. Further, the CPU 210 may communicate witha hardware controller for devices such as for a display 206. The display206, for example, may be used to display text and graphics. One exampleof a suitable display 206 is a touchscreen that provides graphical andtextual visual feedback to a user. In some embodiments, the display 206includes the input devices 204 as part of the display, such as when theinput device is a touchscreen. In some embodiments, the display 206 isseparate from the input device 204. Examples of standalone displaydevices include, for example, an LCD display screen, an LED displayscreen, a projected display (such as a heads-up display device), and soon. Other I/O devices 208 may also be coupled to the CPU 210, such as avideo or audio card, USB or other external devices, printer, speakers,CD-ROM drive, DVD drive, disk drives, Blu-Ray devices, batteryconnection cables, or battery measurement tools. In someimplementations, other I/O devices 208 also include a communicationdevice capable of communicating wirelessly or wire-based with a networknode. The communication device may communicate with another device or aserver through a network using, for example, TCP/IP protocols.

The CPU 210 can access a memory 218. The memory 218 can include one ormore hardware devices for volatile and non-volatile storage, and mayinclude both read-only and writable memory. For example, the memory 218may comprise random access memory (RAM), read-only memory (ROM),writable non-volatile memory, such as flash memory, hard drives, floppydisks, CDs, DVDs, magnetic storage devices, tape drives, device buffers,and so forth. The memory 218 can include non-transitory electricalsignals on the underlying hardware. The memory 218 can include programmemory 212 that contains programs and software, such as an operatingsystem 214 and other application programs 216. The memory 218 alsoincludes data memory 220 that includes any configuration data, settings,user options and preferences that may be needed by the program memory212. The controller 202 may include general purpose or special purposecomputing system environments or configurations. In some embodiments,the controller 202 may not include the illustrated input devices 204 ordisplay 206, but may instead be a component that is accessibleprogrammatically only.

Many embodiments of the technology described below may take the form ofcomputer- or controller-executable instructions, including routinesexecuted by a programmable computer or controller. Those skilled in therelevant art will appreciate that the technology can be practiced oncomputer/controller systems other than those shown and described below.The technology can be embodied in a special-purpose computer, controlleror data processor that is specifically programmed, configured orconstructed to perform one or more of the computer-executableinstructions described below. The technology can also be practiced indistributed environments, where tasks or modules are performed by remoteprocessing devices that are linked through a communications network. Ina distributed computing environment, program modules or subroutines maybe located in local and remote memory storage devices. Aspects of thetechnology described below may be stored or distributed onnon-transitory computer-readable media, including magnetic or opticallyreadable or removable computer disks, as well as distributedelectronically over networks. Data structures and transmissions of dataparticular to aspects of the technology are also encompassed within thescope of the embodiments of the technology.

FIG. 3A and FIG. 3B are block diagrams that illustrate non-limitingexample embodiments of the present disclosure as they may be deployed indevices. FIG. 3A illustrates an example of a self-contained rechargeablebattery device 302. The rechargeable battery device 302 isself-contained in that the controller 110, programmable chip 112, andbattery 114 are all enclosed within a single housing. Such an embodimentwould include a controller 110 that has adequate computing power todetermine the model for use by the programmable chip 112 in controllingthe battery 114 as described fully below.

FIG. 3B illustrates an example of a rechargeable battery device 310 thatis not self-contained in the same way as the rechargeable battery device302)) illustrated in FIG. 3A. As shown, the rechargeable battery device310 has a housing that encloses the battery 114 and the programmablechip 112, but not the controller 110. Instead, the housing of therechargeable battery device 310 encloses a communication interface 306,which is configured to communicate with a corresponding communicationinterface 304 of a model computation system 308. The model computationsystem 308 and controller 110 may be implemented using any suitablecomputing device or devices, including but not limited to a laptopcomputing device, a desktop computing device, a cluster of computingdevices, and one or more computing devices in a cloud computing service.

Communication between the communication interface 306 of therechargeable battery device 310 and the communication interface 304 ofthe model computation system 308 may be by any suitable communicationtechnique, including but not limited to wired or wireless Internetcommunication; wireless communication including but not limited to WiFi,3G, 4G, LTE, or Bluetooth; and wired communication including but notlimited to USB, Firewire, Ethernet, fiber optic, CAN bus, or OBD-II. Inuse, the rechargeable battery device 310 may measure a discharge curvefor the battery 114, and may transmit the discharge curve and any otherrelevant information (including but not limited to temperature readingsor other readings from one or more sensors) to the model computationsystem 308. The model computation system 308 may then use its superiorcomputing power to have the controller 110 determine the model of thebattery 114 based on the discharge curve and any other relevantinformation as described below. The model is then transmitted back tothe rechargeable battery device 310.

FIG. 4 illustrates the computational schematic of the pseudo2-dimensional (p2D) model of Newman and co-workers, a continuumelectrochemical model that has found substantial application forsimulation of Li-ion battery performance. The active material in bothelectrodes is modeled as spherical particles. Electron-transferreactions are modeled at the particle-electrolyte interface, as is thetransport of intercalated lithium through the active particle. Liquidphase mass and charge transport through the thickness of each domainalso modeled using concentrated solution theory. Electron transportthrough the solid phase is also considered, with electronic currententering and leaving the cell at the current collectors (not shown).

The typical p2D model is written for a single ‘cathode-separator-anode’sandwich. Each domain is modeled using porous electrode theory, in whichthe two solid and electrolyte phases are regarded as superimposedcontinua. The model is a set of coupled partial differential equations(PDEs) based on one-dimensional conservation laws for charge and mass ineach domain. The individual domain equations are coupled through thespecification of appropriate interfacial boundary conditions, which alsoensure mathematical well-posedness. In the p2D representation of thebattery, the active material is regarded as composed of sphericalparticles of uniform radii. Lithium intercalation and de-intercalationoccurs through electron-transfer reactions at the particle surface andtransport through the solid particle, modeled by conservation laws inthe ‘pseudo’ r-dimension. The complete mathematical model and parametervalues may be found in Table I-Table III below.

In the absence of complexities such as phase-separation or concentratedsolution effects, solid phase transport is modeled by Fick's second lawin spherical coordinates. For the positive electrode particle, we have

$\begin{matrix}{\frac{\partial c_{1}^{s}}{\partial t} = {\frac{1}{r^{2}}\frac{\partial}{\partial r}{\left( {r^{2}D_{1}^{s}\frac{\partial c_{1}^{s}}{\partial r}} \right)\mspace{31mu}\left\lbrack {0 < x < l_{1}} \right\rbrack}}} & (1)\end{matrix}$

Where the r coordinate denotes the radial distance within the particleand is thus the second ‘pseudo-dimension’. The subscript 1 denotesvariables and parameters pertinent to the positive electrode (region 1).

In porous electrode theory, equation (1) is valid at each point alongthe electrode thickness x. The superscript s is used to denote the solidphase. D₁ ^(s) is the diffusion coefficient in the positive particle.The second-order problem in r requires the specification of two boundaryconditions. At the surface of the solid particle, the diffusive flux isrelated to the local rate of electrode reaction or pore-wall flux as

$\begin{matrix}{{{- D_{1}^{s}}\frac{\partial c_{1}^{s}}{\partial r}} = {j_{1}\mspace{31mu}\left\lbrack {{r = R_{1}},{0 < x < l_{1}}} \right\rbrack}} & (2)\end{matrix}$

A symmetry boundary condition is applied at the center of the positiveparticle

$\begin{matrix}{\frac{\partial c_{1}^{s}}{\partial r} = {0\mspace{31mu}\left\lbrack {{r = 0},{0 < x < l_{1}}} \right\rbrack}} & (3)\end{matrix}$

The analogous set of equations are written for the negative electrode(region 3)

$\begin{matrix}{\frac{\partial c_{3}^{s}}{\partial t} = {\frac{1}{r^{2}}\frac{\partial}{\partial r}{\left( {r^{2}D_{3}^{s}\frac{\partial c_{3}^{s}}{\partial r}} \right)\mspace{31mu}\left\lbrack {{l_{1} + l_{2}} < x < {l_{1} + l_{2} + l_{3}}} \right\rbrack}}} & (4)\end{matrix}$

With the boundary conditions given by

$\begin{matrix}{{{- D_{3}^{s}}\frac{\partial c_{3}^{s}}{\partial r}} = {j_{3}\mspace{31mu}\left\lbrack {{r = R_{3}},{{l_{2} + l_{3}} < x < {l_{1} + l_{2} + l_{3}}}} \right\rbrack}} & {{~~~~~~~~~~~~~~~~~~~~~~~~~~~~~}(5)} \\{\frac{\partial c_{3}^{s}}{\partial r} = {0\mspace{31mu}\left\lbrack {{r = 0},{{l_{1} + l_{2}} < x < {l_{1} + l_{2} + l_{3}}}} \right\rbrack}} & {{~~~~~~~~~~~~~~~~~~~~~~~~~~~~~}(6)}\end{matrix}$

The Butler-Volmer equation is a common constitutive relation for thepore-wall flux in each electrode as follows

$\begin{matrix}{j_{1} = {{\quad\quad}{\quad{{\quad\quad}{\quad{k_{1}{{c_{1}}^{\alpha_{a,1}}\left( {c_{1}^{s,\max} - c_{1}^{s,{surf}}} \right)}^{\alpha_{a,1}}\left( c_{1}^{s,{surf}} \right)^{\alpha_{c,1}}\left( {{\exp\left( \frac{\alpha_{a,1}F\;\eta_{1}}{RT} \right)} - \left. \quad{\exp\left( \frac{{- \alpha_{c,1}}F\;\eta_{1}}{RT} \right)} \right)} \right.}}}}}} & (7)\end{matrix}$

Where k₁ is the rate constant for the positive electrode reaction, andc₁ ^(s,max) denotes the maximum concentration of Li in the positiveelectrode particle. F and R denote Faraday's constant and the gasconstant respectively. α's are the charge transfer coefficients for eachelectrode reaction. The quantity η=ϕ_(s,1)−ϕ_(l,1)−U₁(c₁ ^(s,surf)) isthe surface overpotential, expressed as the difference between the solidand liquid phase potentials minus the electrode open circuit potentialU₁(c₁ ^(s,surf)) vs. a Li/Li⁺ reference electrode. c₁ ^(s,surf) is thesolid particle concentration evaluated at the surface of the particle,i.e.c ₁ ^(s,surf) =c ₁ ^(s)(r=R ₁ ,x,t)  (8)

The dependence of open circuit potential (OCP) on the surfaceconcentration is indicated accordingly. The equivalent expression forthe negative electrode is given by

$\begin{matrix}{{{j_{3} =}\quad}k_{3}{{c_{3}}^{\alpha_{a,3}}\left( {c_{3}^{s,\max} - c_{3}^{s,{surf}}} \right)}^{\alpha_{a,3}}\left( c_{3}^{s,{surf}} \right)^{\alpha_{c,3}}\left( {{\exp\left( \frac{\alpha_{a,3}F\;\eta_{3}}{RT} \right)} - \left. \quad{{\quad\quad}{\exp\left( \frac{{- \alpha_{c,3}}F\;\eta_{3}}{RT} \right)}} \right)} \right.} & (9)\end{matrix}$

Equations (1) and (4) can be volume-averaged over their respectiveelectrode volumes. For the positive electrode, we obtain

$\begin{matrix}{\frac{\partial{\int_{x = 0}^{x = l_{1}}{c_{1}^{s}{dx}}}}{\partial t} = {\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {r^{2}D_{1}^{s}\frac{\partial{\int_{x = 0}^{x = l_{1}}{c_{1}^{s}{dx}}}}{\partial r}} \right)}} & (10)\end{matrix}$

The volume-averaged form becomes

$\begin{matrix}{\frac{\partial\overset{\_}{c_{1}^{s}}}{\partial t} = {\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {r^{2}D_{1}^{s}\frac{\partial\overset{\_}{c_{1}^{s}}}{\partial r}} \right)}} & (11)\end{matrix}$

Similarly, for the negative electrode, we have

$\begin{matrix}{\frac{\partial\overset{\_}{c_{3}^{s}}}{\partial t} = {\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {r^{2}D_{3}^{s}\frac{\partial\overset{\_}{c_{3}^{s}}}{\partial r}} \right)}} & (12)\end{matrix}$

The corresponding boundary conditions can also be expressed involume-averaged form.

$\begin{matrix}{\begin{matrix}{{\frac{\partial\overset{\_}{c_{i}^{s}}}{\partial r} = {- \frac{\overset{\_}{j_{i}}}{D_{i}^{s}}}},{r = R_{i}}} \\{{\frac{\partial\overset{\_}{c_{i}^{s}}}{\partial r} = 0},{r = 0}}\end{matrix}i\;\epsilon\left\{ {1,3} \right\}} & (13)\end{matrix}$

Numerical solution of these equations entails spatial discretization inthe spherical dimension. Discretization of equations (11) and (12)results in a system of Differential Algebraic Equations (DAEs), with aconvenient linear form for constant D_(i) ^(s). For discretization, weemploy an efficient three-parameter model based on a biquadratic profilefor the radial dependence of c_(i) ^(s). This approximation is expectedto ensure higher accuracy than a two-parameter parabolic profile even atrelatively high rates of discharge. The discretized system of equationsis therefore

$\begin{matrix}{\frac{d\overset{\_}{c_{i}^{s,{avg}}}}{dt} = {{- 3}\frac{\overset{\_}{j_{i}}}{R_{i}}}} & (14) \\{\frac{d\overset{\_}{q_{i}^{avg}}}{dt} = {{{- 30}\frac{D_{i}^{s}}{R_{i}^{2}}\overset{\_}{q_{i}^{avg}}} - {\frac{45}{2}\frac{\overset{\_}{j_{i}}}{R_{i}^{2}}}}} & (15) \\{{{{35{\frac{D_{i}^{s}}{R_{i}}\left\lbrack {\overset{\_}{c_{i}^{s,{surf}}} - \overset{\_}{c_{i}^{s,{avg}}}} \right\rbrack}} - {8D_{i}^{s}\overset{\_}{q_{i}^{avg}}}} = {{{- \overset{\_}{j_{i}}}\mspace{34mu} i} \in \left\{ {1,3} \right\}}}\mspace{31mu}} & (16)\end{matrix}$

Where the three-parameter model has been expressed in terms of theparticle-averaged solid phase concentration c_(i) ^(s,avg) , theparticle-averaged concentration gradient q_(i) ^(avg) , and the particlesurface concentration c_(i) ^(s,surf) . The particle averageconcentrations are related to the State of Charge (SoC) at the celllevel and is directly obtained from the simulation results in the aboveformulation.

As mentioned above, this disclosure provides efficient equations for theelectrolyte phase, and therefore the most commonly adopted solid-phaseapproximation is used herein. Importantly, the accuracy of the quarticprofile approximation was quantitatively verified against nearlyerror-free numerical methods (collocation, finite difference) for thecases considered herein.

For charge transport in the solid phase, the governing equation may bewritten as a conservation law for charge as follows:

$\begin{matrix}{{{- \frac{\partial i_{s,1}}{\partial x}} - {a_{1}{Fj}_{1}}} = {0\mspace{31mu}\left\lbrack {0 < x < l_{1}} \right\rbrack}} & (17)\end{matrix}$

The time-derivative for charge density is ignored due toelectroneutrality.

Similarly, we have, for the negative electrode

$\begin{matrix}{{{- \frac{\partial i_{s,3}}{\partial x}} - {a_{3}{Fj}_{3}}} = {0\left\lbrack {{l_{1} + l_{2}} < x < {l_{1} + l_{2} + l_{3}}} \right\rbrack}} & (18)\end{matrix}$

Here, are i_(s,1) and i_(s,3) denote the solid phase current densities.The constitutive equation for the solid phase current density is anOhm's law expression based on the effective electronic conductivity andlocal potential gradient as follows

$\begin{matrix}{{i_{s,1} = {{- {\sigma_{1}\left( {1 - \varepsilon_{1} - \varepsilon_{f,1}} \right)}}\frac{\partial\phi_{s,1}}{\partial x}}}{i_{s,3} = {{- {\sigma_{3}\left( {1 - \varepsilon_{3} - \varepsilon_{f,3}} \right)}}\frac{\partial\phi_{s,3}}{\partial x}}}} & (19)\end{matrix}$

Where 1−ε_(i)−ε_(f,i) is the fraction of solid phase in electrode i,after subtracting the liquid and inert volume fractions. This factorcorrects for the actual conduction pathways in the electrode material.More detailed correction factors may also be applied. The solid phasevolume fraction also appears in the specific interfacial area, which forperfectly spherical particles is given by

$a_{i} = {\frac{3}{R_{i}}{\left( {1 - \varepsilon_{i} - \varepsilon_{fi}} \right).}}$

Volume-averaging gives us

$\begin{matrix}{\frac{i_{1,{x = 0}} - i_{1,{x = l_{1}}}}{l_{1}} = {a_{1}F\overset{\_}{j_{1}}}} & (20)\end{matrix}$

Using the boundary conditions that impose the solid phase currentdensity at the interfaces, one can simplify equation (20) as

$\begin{matrix}{\frac{i_{app}}{l_{1}} = {a_{1}F\overset{\_}{j_{1}}}} & (21)\end{matrix}$

Volume averaging thus connects the applied cell current density i_(app)and average pore-wall flux. The sign convention for the model is soadopted that i_(app) is negative during discharge.

$\begin{matrix}{\overset{\_}{j_{1}} = \frac{i_{app}}{a_{1}{Fl}_{1}}} & (22)\end{matrix}$

Similarly, for the negative electrode, we have

$\begin{matrix}{\overset{\_}{j_{3}} = {- \frac{i_{app}}{a_{3}{Fl}_{3}}}} & (23)\end{matrix}$

Equations (21) and (22) can be used in conjunction with thevolume—averaged forms of equations (1)-(6) to determine the temporalevolution of average solid-phase concentrations for a given i_(app). Ineffect, the active material in each electrode is now modeled by a singlerepresentative particle, the Li concentration profiles through whichwill be influenced by factors such as the applied current densityi_(app), solid phase diffusion coefficients D_(i) ^(s) and thecharacteristic particle radius R_(i). This set of equations alsodetermines the evolution of the averaged surface particle concentration,which in turn affects the surface overpotentials and thus the cellvoltage response through constitutive equations (7) and (9). Theseequations are volume-averaged in order to obtain a relationship betweenthe average pore-wall fluxes j_(i) and the average potentials ϕ_(s,i)and ϕ_(e,i) . This is illustrated for the positive electrode in equation(24) below

$\begin{matrix}{\overset{\_}{j_{1}} = {\frac{\int_{x = 0}^{x = l_{1}}{j_{1}{dx}}}{\int_{x = 0}^{x = l_{1}}{dx}} = \frac{\int_{x = 0}^{x = l_{1}}{k_{1}{c_{1}^{\alpha_{a,1}}\left( {c_{1,\max}^{s} - c_{1}^{s,{surf}}} \right)}^{\alpha_{a,1}}\left( c_{1}^{s,{surf}} \right)^{\alpha_{c,1}}\left( {{\exp\left( \frac{\alpha_{a,1}F\eta_{1}}{RT} \right)} - {\exp\left( \frac{{- \alpha_{c,1}}F\eta_{1}}{RT} \right)}} \right){dx}}}{\int_{x = 0}^{x = l_{1}}{dx}}}} & (24)\end{matrix}$

Unlike the equations for electrolyte concentration, the highlynon-linear nature of the constitutive equation renders evaluation ofequation (24) cumbersome and likely impossible without the use of afull-order solution that gives the actual spatial dependence of thevariables. To this end, averages of non-linear quantities areapproximated by their value at the average values of the variables onwhich they depend (i.e. f(X)≈f(X)).

Mathematically, this can be stated as

$\begin{matrix}{{\overset{\_}{j_{i}}▯{k_{i}\left( \overset{\_}{c_{i}} \right)}^{\alpha_{a,i}}\text{⁠}\left( {c_{i}^{s,\max} - \overset{\_}{c_{i}^{s,{surf}}}} \right)^{\alpha_{a,i}}\text{⁠}\left( \overset{\_}{c_{i}^{s,{surf}}} \right)^{\alpha_{c,i}}\text{⁠}\left( {{\exp\left( \frac{\alpha_{a,i}F\overset{\_}{\eta_{i}}}{RT} \right)} - {\exp\left( \frac{{- \alpha_{c,i}}F\overset{\_}{\eta_{i}}}{RT} \right)}} \right)}{\overset{\_}{\eta_{i}} = {{\overset{\_}{\phi_{1,i}} - \overset{\_}{\phi_{2,i}} - {{U_{i}\left( \overset{\_}{c_{i}^{s,{surf}}} \right)}i}} \in \left\{ {1,3} \right\}}}} & (25)\end{matrix}$

The classic Single Particle Model (SPM) employs an additionalsimplification by ignoring the dynamics of the electrolyte phase.Therefore, c_(i) =c₀, implying that the electrolyte concentration isalways equal to its initial value. Neglecting liquid phase variationsalso means that ϕ_(2,i) is often set to a constant reference, .e.g.ϕ_(2,i) =0 for all i.

Neglect of ohmic and electrolyte concentration effects restricts theaccuracy of SPM to operating regimes characterized by low ohmic losses,low currents, and kinetically limited electrodes, which usually resultin spatially uniform pore-wall flux distributions. To this end, theTanks-in-Series descriptions of electrolyte dynamics are expected toaugment and improve the practical applicability of SPM.

FIG. 5 illustrates a visualization of a non-limiting example embodimentof the Tank Model approach depicting the mass and charge flows in to andout of each ‘tank’ according to various aspects of the presentdisclosure. Since both the electrolyte flux and liquid phase currentdensity are zero at the current collectors, net flows into the liquidphase of the positive ‘tank’ and out of the negative ‘tank’ are zero.The interfacial boundary conditions at the separator define the mass andcharge flows for the middle separator ‘tank’. The electronic currentcarried by the solid phase at the current collectors is denoted bydotted lines. The solid and liquid phases exchange mass and charge at arate determined by the pore-wall flux j_(i) , an internal exchange whichis not shown here. The sign convention is so adopted that i_(app) isnegative during discharge.

In deriving the Tanks-in-Series model, the volume-averaging procedure isapplied first to the solid phase conservation equations, illustratinghow SPM is recovered under certain assumptions. The concept is thenapplied to the electrolyte transport equations, thereby resulting inaveraged equations for the liquid phase ‘Tanks’ of FIG. 5 .

FIG. 6 is a flowchart that illustrates a non-limiting example embodimentof a method of deriving a tanks-in-series model that represents abattery according to various aspects of the present disclosure. Themethod 600 is a high-level description of the general steps performed inderiving the tanks-in-series model. A rigorous description of thespecific steps and computations involved is provided after thediscussion of the general flowchart. Any suitable computing device maybe used to conduct the method 600, including but not limited to thecontroller 110 described above.

From a start block, the method 600 proceeds to block 602, where acomputing device receives a set of parameters that describe the battery114. The parameters of the set of parameters describe physical elementsthat make up the battery 114, and may include one or more of a porosity,an electrode filler fraction, a Bruggeman tortuosity correction, aparticle surface area per unit volume, a maximum particle phaseconcentration, an initial particle phase concentration an initialelectrolyte concentration, a solid phase diffusivity, an electrodereaction rate constant, an electrode reaction anodic coefficient, anelectrode reaction cathodic coefficient, an electrode thickness, acharacteristic particle radius, a Li+ transference number, an electronicconductivity, a temperature, and a current density (1C).

At block 604, the computing device determines a plurality of domains forthe battery 114 based on the set of parameters. In some embodiments,different parameters may be provided for each domain (e.g., an anodedomain, a cathode domain, and a separator domain) of the battery 114,and so the computing device may determine the plurality of domains byreviewing the domains for which parameters were provided. Typically,three domains (an anode domain, a cathode domain, and a separatordomain) are provided for a given battery 114, but this should not beseen as limiting. In some embodiments wherein multiple cells are stackedto form a single battery 114, domains may be provided for each cell.

The method 600 then proceeds to a for-loop defined between a for-loopstart block 606 and a for-loop end block 610 that is executed for eachdomain that was determined at block 604. From the for-loop start block606, the method 600 advances to block 608, where the computing deviceintegrates a pseudo-two dimensional (p2D) model representation for avolume of the domain to create a tank model of the domain. Thisintegration simplifies the calculation that would otherwise be used tosolve the p2D model, thus improving the speed of the calculation andenabling the calculation to be performed on less-capable computinghardware. A detailed description of this integration is provided below.

The method 600 then advances to the for-loop end block 610. If furtherdomains remain to be processed, then the method 600 returns to for-loopstart block 606 to process the next domain. Otherwise, if all of thedomains have been processed, then the method 600 advances to a for-loopdefined between a for-loop start block 612 and a for-loop end block 616that is executed for each interface between domains determined at block604. Typically, an interface will be determined between the anode domainand the separator domain, and another interface will be determinedbetween the separator domain and the cathode domain, though otherarrangements could be possible. From the for-loop start block 612, themethod 600 advances to block 614, where the computing deviceapproximates an interfacial flux for the interface, where thecorresponding gradient is expressed as an average value minus aninterfacial value of a dependent variable over a length scale in amanner that preserves mass and charge conservation between each tankmodel. A detailed description of the calculation and use of theinterfacial flux values is provided below.

The method 600 then advances to the for-loop end block 616. If furtherinterfaces remain to be processed, then the method 600 returns tofor-loop start block 612 to process the next interface. Otherwise, ifall of the interfaces have been processed, then the method 600 advancesto block 618.

At block 618, the computing device combines the tank models and theapproximated interfacial fluxes to create the tanks-in-series model.Because mass and charge have been conserved between each tank model, thetanks-in-series model has a high degree of accuracy becauseconcentration-dependent transport properties may be considered. Further,this combination reduces the full p2D model into a fixed size systemwithout partial differential equations, thus allowing for greatlyincreased speed and greatly reduced processing requirements. A detaileddescription of the combination of the tank models and the approximatedinterfacial fluxes is provided below.

The method 600 then proceeds to an end block and terminates.

The tanks-in-series model computed in the method 600 may be stored in acomputer-readable medium, or may be transmitted by the computing deviceto another computing device for use in controlling charging and/ordischarging of a battery 114. In some embodiments, the tanks-in-seriesmodel may be provided to a controller 110 or a programmable chip 112 inorder to control the charging and/or discharging of a battery 114. FIG.7 is a flowchart that illustrates a non-limiting example embodiment of amethod of controlling at least one of charging and discharging of abattery according to various aspects of the present disclosure. In themethod 700 as illustrated and described, it is assumed that thecontroller 110 has access to the tanks-in-series model computed by themethod 600. In some embodiments, instead of the controller 110, themethod 700 may be performed by the programmable chip 112 after thetanks-in-series model is provided to it by the controller 110.

From a start block, the method 700 proceeds to block 702, where acontroller 110 receives at least one battery state of the battery 114from at least one sensor. The at least one battery state may be anybattery state the tanks-in-series model is designed to accept as input.For example, the at least one battery state may include a terminalvoltage of the battery 114 or a surface temperature of the battery 114sensed by an appropriate sensor device operatively coupled to thebattery 114.

At block 704, the controller 110 provides the at least one battery stateas input to a tanks-in-series model that represents the battery 114. Thetanks-in-series model uses the at least one battery state to determine astate-of-charge, a state of health, or some other value that can be usedto control charging and/or discharging of the battery 114. At block 706,the controller 110 uses at least one output of the tanks-in-series modelto control the at least one of charging and discharging of the battery114.

The method 700 then proceeds to an end block and terminates.

As discussed above, a complete description of how to derive thetanks-in-series model follows.

We begin with the governing equation for electrolyte concentration foran isothermal model in one spatial dimension. The equations for c basedon porous electrode theory may be expressed in the form of conservationlaws

In the positive electrode,

$\begin{matrix}{{\varepsilon_{1}\frac{\partial c_{1}}{\partial t}} = {{- \frac{\partial N_{1}}{\partial x}} + {{a_{1}\left( {1 - t_{+}^{0}} \right)}{j_{1}\left\lbrack {0 < x < l_{1}} \right\rbrack}}}} & (26)\end{matrix}$

Due to the absence of solid active material, the conservation equationfor c in the separator is characterized by a lack of a source term as

$\begin{matrix}{{\varepsilon_{2}\frac{\partial c_{2}}{\partial t}} = {- {\frac{\partial N_{2}}{\partial x}\left\lbrack {l_{1} < x < {l_{1} + l_{2}}} \right\rbrack}}} & (27)\end{matrix}$

The subscript 2 is used to denote the variables in the separator domain.

The equation for the negative electrode is identical in form to that ofthe positive electrode

$\begin{matrix}{{\varepsilon_{3}\frac{\partial c_{3}}{\partial t}} = {{- \frac{\partial N_{3}}{\partial x}} + {{a_{3}\left( {1 - t_{+}^{0}} \right)}{j_{3}\left\lbrack {{l_{1} + l_{2}} < x < {l_{1} + l_{2} + l_{3}}} \right\rbrack}}}} & (28)\end{matrix}$

N₁, N₂, N₃ may regarded as electrolyte fluxes, which need to be relatedto local concentration gradients. Noting the similarity of the governingequations to Fick's second law, we have the following constitutiveequations

$\begin{matrix}{{N_{1} = {{- {D\left( c_{1} \right)}}\varepsilon_{1}^{b_{1}}\frac{\partial c_{1}}{\partial x}}}{N_{2} = {{- {D\left( c_{2} \right)}}\varepsilon_{2}^{b_{2}}\frac{\partial c_{2}}{\partial x}}}{N_{3} = {{- {D\left( c_{3} \right)}}\varepsilon_{3}^{b_{3}}\frac{\partial c_{3}}{\partial x}}}} & (29)\end{matrix}$

In the above equations, D(c) denotes the concentration-dependentelectrolyte diffusion coefficient, corrected by a Bruggeman-type factorto account for porous medium tortuosity.

The governing equations for electrolyte concentration are second-orderin space, which entails the specification of two boundary conditions forc₁, c₂ and c₃. The boundary conditions are defined at the extremities ofeach domain. The positive and negative current collectors are physicalbarriers to the transport of Li⁺ ions, and thus the electrolyte flux atthese locations is set to zero. These boundary conditions are thus givenbyN _(1,x=0) =N ₀₁=0  (30)And,N _(3,x=l) _(p) _(+l) _(s) _(+l) _(n) =N ₃₄=0  (31)

In addition, electrolyte concentrations and their fluxes must becontinuous at the interface between the separator and electrodes. At thepositive electrode-separator interface, this is expressed asc _(1,x=l) ₁ =c _(2,x=l) ₁N _(1,x=l) ₁ =N _(2,x=l) ₁ =N ₁₂  (32)

In general N_(ij) is used to denote the flux at the interface of regionsi and j.

Similarly, at the interface between the separator and negativeelectrode, we havec _(2,x=l) ₁ _(+l) ₂ =c _(3,x=l) ₁ _(+l) ₂N _(2,x=l) ₁ _(+l) ₂ =N _(3,x=l) ₁ _(+l) ₂ =N ₂₃  (33)

Now, equation (26) is integrated over the volume of the positiveelectrode V₁ as

$\begin{matrix}{\frac{\partial{\int_{V_{1}}{\varepsilon_{1}c_{1}{dV}}}}{\partial t} = {{- {\int_{V_{1}}{\frac{\partial N_{1}}{\partial x}{dV}}}} + {\int_{V_{1}}{{a_{1}\left( {1 - t_{+}^{0}} \right)}j_{1}{dV}}}}} & (34)\end{matrix}$

For the one-dimensional model in cartesian coordinates, the differentialvolume dV is given by dV=Adx, where A is a constant that may beconsidered a cross-sectional area. In addition, we express the integralsin terms of average quantities as

${{\overset{\_}{c_{1}} = \frac{\int_{V_{1}}{c_{1}{dV}}}{\int_{V_{1}}{dV}}},{and}}{{\overset{\_}{j_{1}} = \frac{\int_{V_{1}}{j_{1}{dV}}}{\int_{V_{1}}{dV}}},}$with v denoting the volume average of variable v in a given ‘tank’.Substituting these relations converts the volume integral into aone-dimensional integral over electrode thickness, i.e. from x=0 tox=l₁.

Equation (34) thus becomes

$\begin{matrix}{{\varepsilon_{1}\frac{d\overset{\_}{c_{1}}}{dt}} = {{{- \frac{\int_{x = 0}^{x = l_{1}}{\frac{\partial N_{1}}{\partial x}dx}}{l_{1}}} + {{a_{1}\left( {1 - t_{+}^{0}} \right)}\overset{\_}{j_{1}}}} = {\frac{N_{1,{x = 0}} - N_{1,{x = l_{1}}}}{l_{1}} + {{a_{1}\left( {1 - t_{+}^{0}} \right)}\overset{\_}{j_{1}}}}}} & (35)\end{matrix}$

Here, we make the reasonable assumption that the electrode porosity ε₁,specific interfacial area a₁ and Li⁺ transport number t₊ ⁰ in theelectrolyte phase are constant in both space and time.

Using equation (30), equation (35) reduces to

$\begin{matrix}{{\varepsilon_{1}\frac{d\overset{\_}{c_{1}}}{dt}} = {\frac{- N_{{1x} = l_{1}}}{l_{1}} + {{a_{1}\left( {1 - t_{+}^{0}} \right)}\overset{\_}{j_{1}}}}} & (36)\end{matrix}$

The same sequence of operations gives us the volume-averaged equationsfor the separator and negative electrode as follows

$\begin{matrix}{{\varepsilon_{2}\frac{d\overset{\_}{c_{2}}}{dt}} = \frac{N_{2,{x = l_{1}}} - N_{2,{x = {l_{1} + l_{2}}}}}{l_{2}}} & (37)\end{matrix}$ $\begin{matrix}{{\varepsilon_{3}\frac{d\overset{\_}{c_{3}}}{dt}} = {\frac{N_{3,{x = {l_{1} + l_{2}}}}}{l_{3}} + {{a_{3}\left( {1 - t_{+}^{0}} \right)}\overset{\_}{j_{3}}}}} & (38)\end{matrix}$

It is worth noting that the steps applied so far represent a rigorousvolume-averaging of the equations in each porous domain, followed by theuse of the boundary conditions to eliminate interfacial flux terms wherepossible. No approximations have been made up to this point.

In order to track the average concentrations in each ‘tank’, we beginwith the volume-averaged concentration equations (36)-(38). Inspectionof these equations reveals the presence of the unknown interfacial fluxterms that require suitable approximations to achieve closure. In doingso, we can exploit the flux boundary conditions (32) and (33), whichestablish the mass flow coupling between adjacent tanks in series.

A simple approximation for the interfacial diffusive flux is in terms ofa ‘driving force’ Δc and a ‘length scale’ approximation δ_(i,ij) withindomain i for the interface between domains i and j. This is analogous tothe ‘diffusion-length’ approach attempted previously, but we dispensewith assumptions on the spatial profiles for c. The ‘driving force’ Δcis expressed in terms of a difference between the average concentrationand the unknown interfacial concentration Δc₁=c₁ −c_(1,x=l) ₁ , and weuse

$\delta_{1,12} = \frac{l_{1}}{2}$as a first approximation. We therefore have, using the constitutiveequations (29)

$\begin{matrix}{N_{1,{x = l_{1}}} = {{{{- {D\left( c_{1,{x = l_{1}}} \right)}}\varepsilon_{1}^{b_{1}}\frac{\partial c_{1}}{\partial x_{x = l_{1}}}} \approx {{D\left( c_{1,{x = l_{1}}} \right)}{\varepsilon_{1}^{b_{1}}\left( \frac{\Delta c_{1}}{\delta_{1,{12}}} \right)}}} = {{D\left( c_{1,{x = l_{1}}} \right)}{\varepsilon_{1}^{b_{1}}\left( \frac{\overset{\_}{c_{1}} - c_{1,{x = l_{1}}}}{\frac{l_{1}}{2}} \right)}}}} & (39)\end{matrix}$

On the separator side of the interface, we assume

${\delta_{2,{12}} = \frac{l_{2}}{2}},$which we use for the flux approximation

$\begin{matrix}{N_{2,{x = l_{1}}} = {{{{- {D\left( c_{2,{x = l_{1}}} \right)}}\varepsilon_{2}^{b_{2}}\frac{\partial c_{2}}{\partial x_{x = l_{2}}}} \approx {{D\left( c_{2,{x = l_{1}}} \right)}{\varepsilon_{2}^{b_{2}}\left( \frac{\Delta c_{2}}{\delta_{2,12}} \right)}}} = {{D\left( c_{2,{x = l_{1}}} \right)}{\varepsilon_{2}^{b_{2}}\left( \frac{{- \overset{\_}{c_{2}}} + c_{2,{x = l_{1}}}}{\frac{l_{2}}{2}} \right)}}}} & (40)\end{matrix}$

An equivalent interpretation of the above flux approximations is that wehave assumed that the concentration at the midpoint of the porous domainequal to the volume average. This approximation is mathematicallyequivalent to Gaussian integration with one point, accurate to l₁ ².

Substitution of the above approximations into the continuity conditionsof equations (32) gives us

$\begin{matrix}{{{D\left( c_{1,{x = l_{1}}} \right)}{\varepsilon_{1}^{b_{1}}\left( \frac{\overset{\_}{c_{1}} - c_{1,{x = l_{1}}}}{\frac{l_{1}}{2}} \right)}} = {{D\left( c_{2,{x = l_{1}}} \right)}{\varepsilon_{2}^{b_{2}}\left( \frac{{- \overset{\_}{c_{2}}} + c_{2,{x = l_{1}}}}{\frac{l_{2}}{2}} \right)}}} & (41)\end{matrix}$

The interfacial concentration is now expressed in terms of tank averagesas

$\begin{matrix}{c_{12} = {c_{1,{x = l_{1}}} = {c_{2,{x = l_{1}}} = \left( \frac{{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}}\overset{\_}{c_{1}}} + {\frac{\varepsilon_{2}^{b_{2}}}{l_{2}}\overset{\_}{c_{2}}}}{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}} + \frac{\varepsilon_{2}^{b_{2}}}{l_{2}}} \right)}}} & (42)\end{matrix}$

In general, we use the notation v_(ij) to denote the value of a givenvariable v at the interface of domains i and j. An identical sequence ofsteps yields the concentration at the separator-negative interface as

$\begin{matrix}{c_{23} = {c_{2,{x = {l_{1} + l_{2}}}} = {c_{3,{x = {l_{1} + l_{2}}}} = \left( \frac{{\frac{\varepsilon_{3}^{b_{3}}}{l_{3}}\overset{\_}{c_{3}}} + {\frac{\varepsilon_{2}^{b_{2}}}{l_{2}}\overset{\_}{c_{2}}}}{\frac{\varepsilon_{3}^{b_{3}}}{l_{3}} + \frac{\varepsilon_{2}^{b_{2}}}{l_{2}}} \right)}}} & (43)\end{matrix}$

Substituting the values of interfacial concentrations back into the fluxapproximations allows us to express the interfacial fluxes in terms oftank-average variables. Thus, we have

$\begin{matrix}{N_{12} = \frac{{- 2}{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}}} & (44)\end{matrix}$ $\begin{matrix}{N_{23} = \frac{{- 2}{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}} & (45)\end{matrix}$

Equations (44) and (45) can thus be inserted into the volume-averagedforms (36)-(38) to obtain a system of Ordinary Differential Equations(ODEs). It is important to equate the approximations for the total flux,and not just the driving forces. This ensures true mass conservation.

In the three ‘tanks’, after substituting the known values of averagepore-wall fluxes, we have

$\begin{matrix}{\frac{d\overset{\_}{c_{1}}}{dt} = {\frac{\frac{2{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}}}{\varepsilon_{1}l_{1}} + {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{F\varepsilon_{1}l_{1}}}}} & (46)\end{matrix}$ $\begin{matrix}{\frac{d\overset{\_}{c_{2}}}{dt} = \frac{\frac{{- 2}{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} + \frac{2{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}}{\varepsilon_{2}l_{2}}} & (47)\end{matrix}$ $\begin{matrix}{\frac{d\overset{\_}{c_{3}}}{dt} = {\frac{\frac{{- 2}{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}}{\varepsilon_{3}l_{3}} - {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{F\varepsilon_{3}l_{3}}}}} & (48)\end{matrix}$

As specified previously, a slight difference between this approach andothers in literature is that it avoids the assumption of uniformreaction rate (given by the pore-wall flux) to solve the PDEs forconcentration, but instead deals in average quantities. This allows theprediction of average concentration trends even when the constantpore-wall flux assumption is not applicable, with the fluxapproximations as the sole source of error. Inspection of equations(46)-(48) also suggests their decoupling from those for otherelectrochemical variables, indicating that they may be solvedindependently as an ODE system, the solutions of which may be used tocompute other relevant quantities during post-processing. While allmodel equations may be simulated simultaneously throughout thisdisclosure, such a segregated approach may be computationally efficientin real-time control or resource-constrained environments and is enabledby the Tanks-in-Series model. Leaving the model in this form allows forthe incorporation of nonlinear diffusivities.

The governing equation for electrolyte current is related to chargeconservation in the electrolyte phase. Thus, we have

$\begin{matrix}{{{- \frac{\partial i_{l,1}}{\partial x}} = {{- F}a_{1}j_{1}}}{{- \frac{\partial i_{l,2}}{\partial x}} = 0}{{- \frac{\partial i_{l,3}}{\partial x}} = {{- F}a_{3}j_{3}}}} & (49)\end{matrix}$

Where the constitutive equation for electrolyte current is given by amodified Ohm's law based on concentrated solution theory. Thus

$\begin{matrix}{{i_{l,1} = {{{- {\kappa\left( c_{1} \right)}}\varepsilon_{1}^{b_{1}}\frac{\partial\phi_{l,1}}{\partial x}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{1}}} \right){\kappa\left( c_{1} \right)}\varepsilon_{1}^{b_{1}}\frac{1}{c_{1}}\frac{\partial c_{1}}{\partial x}}}}{i_{l,2} = {{{- {\kappa\left( c_{2} \right)}}\varepsilon_{2}^{b_{2}}\frac{\partial\phi_{l,2}}{\partial x}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{2}}} \right){\kappa\left( c_{2} \right)}\varepsilon_{2}^{b_{2}}\frac{1}{c_{2}}\frac{\partial c_{2}}{\partial x}}}}{i_{l,3} = {{{- {\kappa\left( c_{3} \right)}}\varepsilon_{3}^{b_{3}}\frac{\partial\phi_{l,3}}{\partial x}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{3}}} \right){\kappa\left( c_{3} \right)}\varepsilon_{3}^{b_{3}}\frac{1}{c_{3}}\frac{\partial c_{3}}{\partial x}}}}} & (50)\end{matrix}$

The concentration-dependent ionic conductivity κ(c) is corrected by atortuosity factor specific to each region.

Now, volume-averaging equation (49) is redundant, ultimately resultingin equation (21) due to the overall charge balance imposed by porouselectrode theory. However, the interfacial boundary conditions fori_(2,i) provide for the estimation of liquid phase ohmic effects. At theinterface between the electrodes and separator, the entire currenti_(app) is carried by the liquid phase, and the solid-phase currentdensity is zero. Using the constitutive equations (50), we thereforehave, at the positive electrode-separator interface

$\begin{matrix}{{i_{l,1,{x = l_{1}}} = {{{{- {\kappa\left( c_{1,{x = l_{1}}} \right)}}\varepsilon_{1}^{b_{1}}\frac{\partial\phi_{l,1}}{\partial x}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{1,{x = l_{1}}} \right)}{\kappa\left( c_{1,{x = l_{1}}} \right)}\varepsilon_{1}^{b_{1}}\frac{1}{c_{1,{x = l_{1}}}}\frac{\partial c_{1}}{\partial x}}} = i_{app}}}{i_{l,2,{x = l_{1}}} = {{{{- {\kappa\left( c_{2,{x = l_{1}}} \right)}}\varepsilon_{2}^{b_{2}}\frac{\partial\phi_{l,2}}{\partial x}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{2,{x = l_{1}}} \right)}{\kappa\left( c_{2,{x = l_{1}}} \right)}\varepsilon_{2}^{b_{2}}\frac{1}{c_{2,{x = l_{1}}}}\frac{\partial c_{2}}{\partial x}}} = i_{app}}}} & (51)\end{matrix}$

Where we have defined the thermodynamic factor as

${v\left( c_{i} \right)} = {1 + {\frac{{\partial\ln}f}{{\partial\ln}c_{i}}.}}$The tank-averaged equations for the electrolyte potential are nowwritten as

$\begin{matrix}{{i_{l,1,{x = l_{1}}} = {{{{- {\kappa\left( c_{1,{x = l_{1}}} \right)}}{\varepsilon_{1}^{b_{1}}\left( \frac{\phi_{l,{x = l_{1}}} - \overset{\_}{\phi_{l,1}}}{\frac{l_{p}}{2}} \right)}} + {\frac{2{{RT}\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{1,{x = l_{1}}} \right)}{\kappa\left( c_{1,{x = l_{1}}} \right)}\varepsilon_{1}^{b_{1}}\frac{1}{c_{1,{x = l_{1}}}}\left( \frac{c_{1,{x = l_{1}}} - \overset{\_}{c_{1}}}{\frac{l_{p}}{2}} \right)}} = i_{app}}}{i_{l,2,{x = l_{1}}} = {{{{- {\kappa\left( c_{2,{x = l_{1}}} \right)}}{\varepsilon_{2}^{b_{2}}\left( \frac{\overset{\_}{\phi_{l,2}} - \phi_{l,{x = l_{1}}}}{\frac{l_{s}}{2}} \right)}} + {\frac{2{{RT}\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{2,{x = l_{1}}} \right)}{\kappa\left( c_{2,{x = l_{1}}} \right)}\varepsilon_{2}^{b_{2}}\frac{1}{c_{2,{x = l_{1}}}}\left( \frac{\overset{\_}{c_{2}} - c_{2,{x = l_{1}}}}{\frac{l_{s}}{2}} \right)}} = i_{app}}}} & (52)\end{matrix}$

Where an approximation

${{\frac{\partial\phi_{l,1}}{\partial x} \approx {\frac{\left( {\phi_{l,1,{x = l_{p}}} - \overset{\_}{\phi_{l,1}}} \right)}{\frac{l_{p}}{2}}{and}\frac{\partial\phi_{l,2}}{\partial x}}} = \frac{\left( {\phi_{l,2,{x = l_{p}}} - \overset{\_}{\phi_{l,2}}} \right)}{\frac{l_{s}}{2}}},$analogous to equations (39) and (40) is made for the gradients ofϕ_(l,i), in terms of tank-average and interfacial values. Equating theinterfacial current density results in equation (53), which can besolved in conjunction with continuity conditions to obtain theinterfacial electrolyte potential.

$\begin{matrix}{{- {\varepsilon_{1}^{b_{1}}\left( \frac{\phi_{l,1,{x = l_{1}}} - \overset{\_}{\phi_{l,1}}}{\frac{l_{1}}{2}} \right)}} = {- {\varepsilon_{2}^{b_{2}}\left( \frac{\overset{\_}{\phi_{l,2}} - \phi_{l,2,{x = l_{1}}}}{\frac{l_{2}}{2}} \right)}}} & (53)\end{matrix}$

Solving for ϕ_(l,12) gives us

$\begin{matrix}{\phi_{l,12} = {\phi_{l,1,{x = l_{1}}} = {\phi_{l,2,{x = l_{1}}} = \left( \frac{{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}}\overset{\_}{\phi_{l,1}}} + {\frac{\varepsilon_{2}^{b_{2}}}{l_{2}}\overset{\_}{\phi_{l,2}}}}{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}} + \frac{\varepsilon_{2}^{b_{2}}}{l_{2}}} \right)}}} & (54)\end{matrix}$

A similar form is obtained at the separator-negative electrode interface

$\begin{matrix}{\phi_{l,23} = {\phi_{l,2,{x = {l_{1} + l_{2}}}} = {\phi_{l,3,{x = {l_{1} + l_{2}}}} = \left( \frac{{\frac{\varepsilon_{3}^{b_{3}}}{l_{3}}\overset{\_}{\phi_{l,3}}} + {\frac{\varepsilon_{2}^{b_{2}}}{l_{2}}\overset{\_}{\phi_{l,2}}}}{\frac{\varepsilon_{3}^{b_{3}}}{l_{3}} + \frac{\varepsilon_{2}^{b_{2}}}{l_{2}}} \right)}}} & (55)\end{matrix}$

Equations (39)-(43) can now be combined with equations (51)-(55) to giveus the algebraic equations governing electrolyte potential

$\begin{matrix}{{i_{app} = {{{{- 2}{\kappa\left( c_{12} \right)}\left( \frac{\overset{\_}{\phi_{l,2}} - \overset{\_}{\phi_{l,1}}}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} \right)} + {\frac{4R{T\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{12} \right)}{\kappa\left( c_{12} \right)}\frac{1}{c_{12}}\left( \frac{\overset{\_}{c_{2}} - \overset{\_}{c_{1}}}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} \right)}} = i_{l,1,{x = l_{1}}}}}{i_{app} = {{{{- 2}{\kappa\left( c_{23} \right)}\left( \frac{\overset{\_}{\phi_{l,3}} - \overset{\_}{\phi_{l,2}}}{\frac{l_{3}}{\varepsilon_{3}^{b_{3}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} \right)} + {\frac{4R{T\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{23} \right)}{\kappa\left( c_{23} \right)}\frac{1}{c_{23}}\left( \frac{\overset{\_}{c_{3}} - \overset{\_}{c_{2}}}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}} \right)}} = i_{l,3,{x = {l_{1} + l_{2}}}}}}} & (56)\end{matrix}$

The final step in the model formulation is the specification of areference potential. A convenient reference is the electrolyte potentialat the interface between the separator and positive electrode. We thusset ϕ_(l,12)=0. This modifies equation (54), and completes the DAEsystem for the Tank Model.

$\begin{matrix}{\phi_{l,12} = {\left( \frac{{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}}\overset{\_}{\phi_{l,1}}} + {\frac{\varepsilon_{2}^{b_{2}}}{l_{2}}\overset{\_}{\phi_{l,2}}}}{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}} + \frac{\varepsilon_{2}^{b_{2}}}{l_{2}}} \right) = 0}} & (57)\end{matrix}$

We can now assemble the complete Tanks-in-Series Model listed below inTable IV.

Solving the system of Table IV now allows the prediction of the cellvoltage asV _(cell)=ϕ_(s,1) −ϕ_(s,3)   (58)

In using equation (58) to calculate cell voltage for the Tanks-in-SeriesModel, it is implicitly assumed that the solid phase potentials at thecell termini may be approximated by their respective electrode averages.In reality, the terminal potentials are determined via interfacialapproximations for the constitutive equations (19) analogous toequations (52) for ϕ_(l,i). Accounting for this potential drop may benecessitated at high current densities, or for electrodes with poorelectronic conductivity. In practice, σ_(eff) ˜50 S/m and the ϕ_(s,i)gradients within the electrode are negligible. This can be seen by thefollowing rough estimation for an aggressive |i_(app)|=100 A/m²

$\begin{matrix}{{{i_{app}▯} - {\sigma_{eff}\frac{\overset{\_}{\phi_{s,1}} - \phi_{s,01}}{l_{1}/2}}}{{\overset{\_}{\phi_{s,1}} - \phi_{s,01}} = {{❘{{- l_{1}}{i_{app}/\left( {2\sigma_{eff}} \right)}}❘} = {{10^{- 2}/2} \sim {{0.0}05V}}}}} & (59)\end{matrix}$

Thus, the upper limit of this solid phase ohmic drop is ˜5 mV perelectrode. This justifies the assumption of uniform ϕ_(s,i) for mostsituations of practical salience.

Table III lists the parameter values for a 1.78 AhNickel-Cobalt-Manganese (NCM)/graphite power cell, collated from varioussources by Tanim et al. Electrolyte transport property correlations weretaken from the work of Valøen and Reimers, and are listed in Table II. Amodified value of R_(i)=1 μm was used for the particle radii, in orderto ensure rapid diffusion with negligible gradients. The absence ofsolid phase diffusion limitations serves the practical purpose ofensuring the accuracy of our three-parameter model for solid phasetransport even at relatively high discharge rates, based on thequantitative guidelines of Subramanian et al. This prevents numericalerrors for solid phase discretization from confounding our analysis,which is focused on examining the accuracy of our Tanks-in-Seriesequations for the liquid phase.

The tanks-in-series may also be recast. Recalling the volume-averagedelectrolyte equations (46)-(48) of the Tank Model, we have

$\begin{matrix}{\frac{d\overset{\_}{c_{1}}}{dt} = {\frac{\frac{2{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}}}{\varepsilon_{1}l_{1}} + {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{F\varepsilon_{1}l_{1}}}}} & (60)\end{matrix}$ $\begin{matrix}{\frac{d\overset{\_}{c_{2}}}{dt} = \frac{\frac{{- 2}{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} + \frac{2{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}}{\varepsilon_{2}l_{2}}} & (61)\end{matrix}$ $\begin{matrix}{\frac{d\overset{\_}{c_{3}}}{dt} = {\frac{\frac{{- 2}{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}}{\varepsilon_{3}l_{3}} - {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{F\varepsilon_{3}l_{3}}}}} & (62)\end{matrix}$

The above equations may be conveniently written as

$\begin{matrix}{\frac{d\overset{\_}{c_{1}}}{dt} = {{\frac{K_{12}}{V_{l,1}}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)} + {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{{FV}_{l,1}}}}} & (63)\end{matrix}$ $\begin{matrix}{\frac{d\overset{\_}{c_{2}}}{dt} = \frac{{- {K_{12}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}} + {K_{23}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}}{V_{l,2}}} & (64)\end{matrix}$ $\begin{matrix}{\frac{d\overset{\_}{c_{3}}}{dt} = {{\frac{- K_{23}}{V_{l,3}}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)} - {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{{FV}_{l,3}}}}} & (65)\end{matrix}$

Where the governing equations have been expressed in terms of equivalent‘mass transfer’ coefficients, volume and source terms. In the aboveequations, the

$K_{ij} = \frac{2{D\left( c_{ij} \right)}}{\frac{l_{i}}{\varepsilon_{i}^{b_{i}}} + \frac{l_{j}}{\varepsilon_{j}^{b_{j}}}}$have units of mass transfer coefficient (m/s), and V_(l,i) denotes theelectrolyte volume (per unit area) in region i with units of length. Thetransport parameters in the electrolyte current balance (56) may also begrouped similarly

$\begin{matrix}{{i_{app} = {{- {S_{12}\left( {\overset{\_}{\phi_{l,2}} - \overset{\_}{\phi_{l,1}}} \right)}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}{S_{12}\left( \frac{\overset{\_}{c_{2}} - \overset{\_}{c_{1}}}{c_{12}} \right)}}}}{i_{app} = {{- {S_{23}\left( {\overset{\_}{\phi_{l,3}} - \overset{\_}{\phi_{l,2}}} \right)}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}{S_{23}\left( \frac{\overset{\_}{c_{3}} - \overset{\_}{c_{2}}}{c_{23}} \right)}}}}} & (66)\end{matrix}$

Where, the thermodynamic factor v(c_(ij))=1 for simplicity. Here,

$S_{ij} = \frac{2{\kappa\left( c_{ij} \right)}}{\frac{l_{i}}{\varepsilon_{i}^{b_{i}}} + \frac{l_{j}}{\varepsilon_{j}^{b_{j}}}}$have units of S/m². The lengths used in the flux approximation, whichfor this case is

$\frac{l_{i}}{2},$are present in the expressions for K_(ij) and S_(ij). Estimating thesequantities is equivalent to deriving a suitable flux or gradientapproximation. For constant diffusivity, the electrolyte concentrationequations constitute a coupled linear system of Ordinary DifferentialEquations with constant coefficients. The form of these equationssuggests a possible parameter estimation strategy in which the TankModel is matched against experimental data by direct estimation of the‘consolidated’ transfer coefficients. Variations due toconcentration-dependent transport properties, or transients in which the‘steady state layer’ is attained over a characteristic timescalecomparable to the total discharge time can potentially be handled byestimating average values over the course of the process.

The full p2D model is used as the benchmark in evaluating thepredictions of the Tanks-in-Series model. The p2D model was discretizedand solved using coordinate transformation, model reformulation andorthogonal collocation techniques described in previous work. The numberof collocation points in each region were adjusted to achieve numericalconvergence of discharge curves at different current densities,ultimately selecting n=(7,3,7) Gauss-Legendre collocation points.Comparisons of the Tank Model with SPM are also reported, for which aFinite-Difference discretization with n=256 internal node points wereemployed. All DAE systems were consistently initialized, and solvedusing the dsolve function in Maple 2018. An absolute solver toleranceabserr=10⁻⁸ was specified.

For evaluating computational performance, the discretized equations weresolved in time using IDA, an Implicit Differential-Algebraic solver inANSI-standard C language under BSD license. IDA is an efficient solverfor initial value problems (IVP) for systems of DAEs, which is part ofthe SUNDIALS (SUite of Nonlinear and DIfferential/ALgebraic equationSolvers) package. The absolute solver tolerance was set to atol=10⁻⁸ anda relative tolerance of rtol=10⁻⁷ was specified.

All computations were performed on an Intel® Core™ i7-7700K processorwith a clock speed of 4.2 GHz, 8 logical cores and 64 GB RAM.

The Tanks-in-Series model was simulated for galvanostatic discharge atC-rates of 1C, 2C, and 5C. In comparing the model predictions to thefull p2D model, we seek to ascertain the accuracy of the Tanks-in-Seriesapproximations, and to quantify its deviation from p2D as a function ofapplied current density. The results were also compared against curvesfrom SPM, chiefly to illustrate improvement over commonly used fastphysics-based models.

The cell voltage-time predictions of our Tanks-in-Series model(hereafter termed the ‘Tank Model’) were compared with SPM and the fullp2D model. The agreement of all three models at 1C and 2C indicatedrelatively low liquid phase polarizations compared to otheroverpotentials, possibly due to the values of the Base Case parameters,which are those of a power-cell rated for higher discharge rates. Evenso, there is a discernible improvement with the Tank Model. The improvedaccuracy of the tank model at 5C discharge is substantially morepronounced, particularly during the intermediate phase of the dischargeprocess. The qualitative trends are identical for both the Tank Modeland SPM, but the curves for the Tank Model appear shifted by nearly aconstant value compared to SPM. This difference is due to the TankModel's estimates of liquid phase ohmic drops and concentrationoverpotentials, which reduces error. The magnitude of theseoverpotentials increases with current density as characterized by theC-rate.

The increased accuracy of the Tank Model at 5C was proven in thesetests, with the absolute error only slightly exceeding 20 mV, incontrast to >60 mV for SPM. A more than threefold reduction in RMSE isobtained at 5C due to the incorporation of the lumped equations for theelectrolyte.

Despite the chemistry-specific considerations and limitations of theTank Model equations, the Tank Model results in an RMSE of 14.3 mV evenat 5C discharge rate, making it competitive in terms of error metricsfor online applications. This error is expected to reduce even furtherbased on the specific battery chemistries being considered, such as inthe case of a negative electrode with slower kinetics and a moremonotonic U(c_(i) ^(s,surf)) curve, which is expected to result in aninherently more uniform reaction profile. This would further bolster thecase for the Tank Model with its substantially improved computationalefficiency.

The Tank Model was compared against SPM because of its substantialubiquity in advanced BMS applications, as SPM has been the model ofchoice when physical detail is desired. In introducing the Tank Model,we seek to propose an alternative for various applications where SPM issubstantially common. However, for additional perspective, the errormetrics of the Tank Model were also considered relative to a recentversion of the SPMe, which corrects SPM with electrolyte dynamics asmentioned previously. The representative values of electrolyte transportproperties required by the SPMe formulation were evaluated at theinitial electrolyte concentration.

The error metrics up to 3C discharge rate were determined. Thesubstantial improvement in error in going from SPM to lumped electrolytemodels was evident, particularly at 3C, where the RMSE reduction isnearly three-fold for SPMe. However, the Tank Model exhibits even lowererror (13 mV vs. 35 mV). This suggests more accurate estimation ofconcentration overpotentials and liquid phase ohmic drops for the chosenparameters. The SPMe uses representative values, whereas the Tank Modelonly utilizes the evaluation of electrolyte transport properties at theinterfaces. Consequently, the concentration-dependent transportproperties can be included in an efficient manner. The expressions inSPMe may need to be rederived to account for these variations.

Computational times for the Tank Model are listed in Table V. For theTank Model, each conservation law in the electrolyte phase is replacedby its volume-averaged form, while the solid phase in each electrode isreplaced by 3 linear DAEs. The original ˜200-1000 DAEs are thus replacedby 14 average conservation equations. This results in a 1C dischargecurve being simulated in ˜2 ms, in contrast to >1000 ms for a standardFinite Difference implementation. The computational speed of this modelis comparable to fast SPM implementations. The Tank Model is alsocompetitive with the state-of-the-art reformulated p2D model from ourgroup, achieving up to an order of magnitude reduction in computationtime. The mathematical similarity of the Tank Model to the reformulatedmodel with n=(1,1,1) collocation prompts the question as to why the(1,1,1) reformulated model is not used instead of developing the TankModel. This is because the proposed model is conservative and exhibitedhigher accuracy. In addition, it can be rewritten so as to contain feweradjustable parameters compared to reformulated models. Increasing theaccuracy and numerical convergence of reformulated models uses itsability to guarantee convergence for different chemistries, parametersand operating conditions by increasing the number of collocation points.

The particulars shown herein are by way of example and for purposes ofillustrative discussion of the preferred embodiments of the presentinvention only and are presented in the cause of providing what isbelieved to be the most useful and readily understood description of theprinciples and conceptual aspects of various embodiments of theinvention. In this regard, no attempt is made to show structural detailsof the invention in more detail than is necessary for the fundamentalunderstanding of the invention, the description taken with the drawingsand/or examples making apparent to those skilled in the art how theseveral forms of the invention may be embodied in practice.

As used herein and unless otherwise indicated, the terms “a” and “an”are taken to mean “one”, “at least one” or “one or more”. Unlessotherwise required by context, singular terms used herein shall includepluralities and plural terms shall include the singular.

Unless the context clearly requires otherwise, throughout thedescription and the claims, the words ‘comprise’, ‘comprising’, and thelike are to be construed in an inclusive sense as opposed to anexclusive or exhaustive sense; that is to say, in the sense of“including, but not limited to”. Words using the singular or pluralnumber also include the plural and singular number, respectively.Additionally, the words “herein,” “above,” and “below” and words ofsimilar import, when used in this application, shall refer to thisapplication as a whole and not to any particular portions of theapplication.

The description of embodiments of the disclosure is not intended to beexhaustive or to limit the disclosure to the precise form disclosed.While the specific embodiments of, and examples for, the disclosure aredescribed herein for illustrative purposes, various equivalentmodifications are possible within the scope of the disclosure, as thoseskilled in the relevant art will recognize.

All of the references cited herein are incorporated by reference.Aspects of the disclosure can be modified, if necessary, to employ thesystems, functions, and concepts of the above references and applicationto provide yet further embodiments of the disclosure. These and otherchanges can be made to the disclosure in light of the detaileddescription.

Specific elements of any foregoing embodiments can be combined orsubstituted for elements in other embodiments. Moreover, the inclusionof specific elements in at least some of these embodiments may beoptional, wherein further embodiments may include one or moreembodiments that specifically exclude one or more of these specificelements. Furthermore, while advantages associated with certainembodiments of the disclosure have been described in the context ofthese embodiments, other embodiments may also exhibit such advantages,and not all embodiments need necessarily exhibit such advantages to fallwithin the scope of the disclosure.

TABLE I Governing PDEs for the p2D model. Equations Boundary ConditionsPositive Electrode (Region 1)${\varepsilon_{1}\frac{\partial c_{1}}{\partial t}} = {{\frac{\partial}{\partial x}\left\lbrack {{D\left( c_{1} \right)}\varepsilon_{1}^{b_{1}}\frac{\partial c_{1}}{\partial x}} \right\rbrack} + {{a_{1}\left( {1 - t_{+}} \right)}j_{1}}}$$\begin{matrix}{\left. \frac{\partial c_{1}}{\partial x} \right|_{x = 0} = 0} \\{{{\varepsilon_{1}^{b_{1}}\frac{\partial c_{1}}{\partial x}}❘_{x = l_{1}}} = {{\varepsilon_{2}^{b_{2}}\frac{\partial c_{2}}{\partial x}}❘_{x = l_{1}}}}\end{matrix}$$i_{l,1} = {{{- {\kappa\left( c_{1} \right)}}\varepsilon_{1}^{b_{1}}\frac{\partial\phi_{l,1}}{\partial x}} + {\frac{2{{RT}\left( {1 - t_{+}^{0}} \right)}}{F}\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{1}}} \right){\kappa\left( c_{1} \right)}\varepsilon_{1}^{b_{1}}\frac{1}{c_{1}}\frac{\partial c_{1}}{\partial x}}}$$\begin{matrix}{\left. \frac{\partial\phi_{l,1}}{\partial x} \right|_{1,{x = 0}} = 0} \\{{{\varepsilon_{1}^{b_{1}}\frac{\partial\phi_{l,1}}{\partial x}}❘_{x = l_{1}}} = {{\varepsilon_{2}^{b_{2}}\frac{\partial\phi_{l,2}}{\partial x}}❘_{x = l_{1}}}}\end{matrix}$${\frac{\partial}{\partial x}\left\lbrack {\sigma_{{eff},1}\frac{\partial\phi_{s,1}}{\partial x}} \right\rbrack} = {a_{1}{Fj}_{1}}$$\begin{matrix}{\left. \frac{\partial\phi_{l,1}}{\partial x} \right|_{x = 0} = {- \frac{i_{app}}{\sigma_{{eff},1}}}} \\{{\frac{\partial\phi_{l,1}}{\partial x}❘_{x = l_{1}}} = 0}\end{matrix}$$\frac{\partial c_{1}^{s}}{\partial t} = {\frac{1}{r^{2}}{\frac{\partial}{\partial r}\left\lbrack {r^{2}D_{1}^{s}\frac{\partial c_{1}^{s}}{\partial x}} \right\rbrack}}$$\begin{matrix}{{\frac{\partial c_{1}^{s}}{\partial r}❘_{r = 0}} = 0} \\{{\frac{\partial c_{1}^{s}}{\partial r}❘_{r = R_{1}}} = {- \frac{j_{1}}{D_{1}^{s}}}}\end{matrix}$ Separator (Region 2)${\varepsilon_{2}\frac{\partial c_{2}}{\partial t}} = {\frac{\partial}{\partial x}\left\lbrack {{D\left( c_{2} \right)}\frac{\partial c_{2}}{\partial x}} \right\rbrack}$$\begin{matrix}{{c_{1}❘_{x = l_{1}}} = {c_{2}❘_{x = l_{1}}}} \\{{c_{2}❘_{x = {l_{1} + l_{2}}}} = {c_{3}❘_{x = {l_{1} + l_{2}}}}}\end{matrix}$$i_{l,2} = {{{- {\kappa\left( c_{2} \right)}}\varepsilon_{2}^{b_{2}}\frac{\partial\phi_{l,2}}{\partial x}} + {\frac{2R{T\left( {1 - t_{+}^{0}} \right)}}{F}\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{2}}} \right){\kappa\left( c_{2} \right)}\varepsilon_{2}^{b_{2}}\frac{1}{c_{2}}\frac{\partial c_{2}}{\partial x}}}$$\begin{matrix}{{\phi_{l,1}❘_{x = l_{1}}} = {\phi_{l,2}❘_{x = l_{1}}}} \\{{\phi_{l,2}❘_{x = {l_{1} + l_{2}}}} = {\phi_{l,3}❘_{x = {l_{1} + l_{2}}}}}\end{matrix}$ Negative Electrode (Region 3)${\varepsilon_{3}\frac{\partial c_{3}}{\partial t}} = {{\frac{\partial}{\partial x}\left\lbrack {{D\left( c_{3} \right)}\frac{\partial c_{3}}{\partial x}} \right\rbrack} + {{a_{3}\left( {1 - t_{+}} \right)}j_{3}}}$$\begin{matrix}{\left. \frac{\partial c_{3}}{\partial x} \right|_{x = {l_{1} + l_{2} + l_{3}}} = 0} \\{\left. {\varepsilon_{2}^{b_{2}}\frac{\partial c_{2}}{\partial x}} \right|_{x = {l_{1} + l_{2}}} = \left. {\varepsilon_{3}^{b_{3}}\frac{\partial c_{3}}{\partial x}} \right|_{x = {l_{1} + l_{2}}}}\end{matrix}$$i_{l,3} = {{{- {\kappa\left( c_{3} \right)}}\varepsilon_{3}^{b_{3}}\frac{\partial\phi_{l,3}}{\partial x}} + {\frac{2{{RT}\left( {1 - t_{+}^{0}} \right)}}{F}\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{3}}} \right){\kappa\left( c_{3} \right)}\varepsilon_{3}^{b_{3}}\frac{1}{c_{3}}\frac{\partial c_{3}}{\partial x}}}$$\begin{matrix}\left. \phi_{l,3} \middle| {}_{x = {l_{1} + l_{2} + l_{3}}}{= 0} \right. \\{\left. {\varepsilon_{2}^{b_{2}}\frac{\partial\phi_{l,2}}{\partial x}} \right|_{x = {l_{1} + l_{2}}} = \left. {\varepsilon_{3}^{b_{3}}\frac{\partial\phi_{l,3}}{\partial x}} \right|_{x = {l_{1} + l_{2}}}}\end{matrix}$${\frac{\partial}{\partial x}\left\lbrack {\sigma_{{eff},3}\frac{\partial\phi_{s,3}}{\partial x}} \right\rbrack} = {a_{3}{Fj}_{3}}$$\begin{matrix}{\left. \frac{\partial\phi_{s,3}}{\partial x} \right|_{x = {l_{1} + l_{2}}} = 0} \\{\left. \frac{\partial\phi_{s,3}}{\partial x} \right|_{x = {l_{1} + l_{2} + l_{3}}} = {- \frac{i_{app}}{\sigma_{{eff},3}}}}\end{matrix}$$\frac{\partial c_{3}^{s}}{\partial t} = {\frac{1}{r^{2}}{\frac{\partial}{\partial r}\left\lbrack {r^{2}D_{3}^{s}\frac{\partial c_{3}^{s}}{\partial r}} \right\rbrack}}$$\begin{matrix}{\left. \frac{\partial c_{3}^{s}}{\partial r} \right|_{r = 0} = 0} \\{\left. \frac{\partial c_{3}^{s}}{\partial r} \right|_{r = R_{3}} = {- \frac{j_{3}}{D_{3}^{s}}}}\end{matrix}$

TABLE II Additional constitutive equations for the p2D model.$j_{1} = {k_{1}{c_{1}^{\alpha_{a,1}}\left( {c_{1}^{s,\max} - c_{1}^{s,{surf}}} \right)}^{\alpha_{a,1}}\left( c_{1}^{x,{surf}} \right)^{\alpha_{c,1}}\left( {{\exp\left( \frac{\alpha_{a,1}F\eta_{1}}{RT} \right)} - {\exp\left( \frac{{- \alpha_{c,1}}F\eta_{1}}{RT} \right)}} \right)}$$j_{3} = {k_{3}{c_{3}^{\alpha_{a,3}}\left( {c_{3}^{s,\max} - c_{3}^{s,{surf}}} \right)}^{\alpha_{a,3}}\left( c_{3}^{s,{surf}} \right)^{\alpha_{c,3}}\left( {{\exp\left( \frac{\alpha_{a,3}F\eta_{3}}{RT} \right)} - {\exp\left( \frac{{- \alpha_{c,3}}F\eta_{3}}{RT} \right)}} \right)}$${{\kappa\left( c_{i} \right)} = {1 \times 10^{- 4}{c_{i}\begin{bmatrix}{\left( {{- 10.5} + {0.0740T} - {{6.9}6 \times 10^{- 5}T^{2}}} \right) +} \\{{c_{i}\left( {{{0.6}68} - {0.0178T} + {2\text{.8} \times 10^{- 5}T^{2}}} \right)} +} \\{c_{i}^{2}\left( {{{0.4}94} - {{8.8}6 \times 10^{- 4}T^{2}}} \right)}\end{bmatrix}}^{2}}},{i \in \left\{ {1,2,3} \right\}}$ σ_(eff,i) = σ_(i)(1 − ε_(i) − ε_(f,i)), i ∈ {1, 3}${{D\left( c_{i} \right)} = {{0.0}001 \times 10^{{- {\lbrack{4.43 + \frac{54}{T - {229} - {0.005c_{i}}}}\rbrack}} - {0.22c_{i}}}}},{i \in \left\{ {1,2,3} \right\}}$${a_{i} = {\frac{3}{R_{i}}\left( {1 - \varepsilon_{i}\  - \varepsilon_{f,i}} \right)}},{i = 1},2,3$U_(p) = −10.72θ_(p) ⁴ + 23.88θ_(p) ³ − 16.77θ_(p) ² + 2.595θ_(p) + 4.563$\theta_{p} = \frac{c_{1}^{s,{surf}}}{c_{1}^{s,\max}}$ U_(n) = 0.1493 +0.8493e^(−61.79θ) ^(n) + 0.3824e^(−665.8θ) ^(n) − e^(39.42θ) ^(n)^(−41.92) − 0.03131 tan⁻¹ (25.59θ_(n) − 4.099) − 0.009434 tan⁻¹(32.49θ_(n) − 15.74)$\theta_{p} = \frac{c_{3}^{s,{surf}}}{c_{3}^{s,\max}}$${{\left( {1 - t_{+}^{0}} \right)\left( {1 + \frac{{\partial\ln}f}{{\partial\ln}c_{i}}} \right)} = {0.601 - {7.5894 \times 10^{- 3}c_{i}^{0.5}} + {3.1053 \times 10^{- 5}\left( {2.5236 - {0.0052T}} \right)c_{i}^{1.5}}}},{i \in \left\{ {1,2,3} \right\}}$

TABLE III Base case model parameters (constant values) Negative PositiveSymbol Parameter Electrode Separator Electrode Units ε_(i) Porosity 0.30.4 0.3 ε_(fi) Electrode filler fraction 0.038 0.12 b_(i) Bruggemantortuosity 1.5 1.5 1.5 correction a_(i) Particle surface area 17400001986000 m²/m³ per unit volume c_(i) ^(s, max) Maximum particle 3108051830 mol/m³ phase concentration c_(i) ^(s, 0) Initial particle phase24578 18645 mol/m³ concentration c₀ Initial electrolyte 1200    mol/m³concentration D_(i) ^(s) Solid phase diffusivity  1.4 × 10⁻¹⁴  2.0 ×10⁻¹⁴ m²/s k_(i) Electrode reaction rate 6.626 × 10⁻¹⁰ 2.405 × 10⁻¹⁰m^(2.5)/(mol^(0.5)s) constant α_(a, i) Electrode reaction 0.5 0.5 anodiccoefficient α_(c, i) Electrode reaction 0.5 0.5 cathodic coefficientl_(i) Electrode thickness   40 × 10⁻⁶  25 × 10⁻⁶ 36.55 × 10⁻⁶ m R_(i)Characteristic particle 10⁻⁶ 10⁻⁶ m radius t₊ ⁰ Li+ transference    0.38number σ_(i) Electronic conductivity 100 100 S/m T Temperature   298.15 K i_(app) Current Density (1 C) −17.54 A/m²

TABLE IV Governing Equations of the Tanks-in-Series Model PositiveElectrode (Region 1)$\frac{d\overset{\_}{c_{1}}}{dt} = {\frac{\frac{2{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}}}{\varepsilon_{1}l_{1}} + {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{F\varepsilon_{1}l_{1}}}}$$i_{app} = {{{- 2}{\kappa\left( c_{12} \right)}\left( \frac{\overset{\_}{\phi_{l,2}} - \overset{\_}{\phi_{l,1}}}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} \right)} + {\frac{4R{T\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{12} \right)}{\kappa\left( c_{12} \right)}\frac{1}{c_{12}}\left( \frac{\overset{\_}{c_{2}} - \overset{\_}{c_{1}}}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} \right)}}$$\frac{d\overset{\_}{c_{1}^{s,{avg}}}}{dt} = {{- 3}\frac{\overset{\_}{j_{1}}}{R_{1}}}$$\frac{d\overset{\_}{q_{1}^{avg}}}{dt} = {{{- 3}0\frac{D_{1}^{s}}{R_{1}^{2}}\overset{\_}{q_{1}^{avg}}} - {\frac{45}{2}\frac{\overset{\_}{j_{1}}}{R_{1}^{2}}}}$${{35{\frac{D_{1}^{s}}{R_{1}}\left\lbrack {\overset{\_}{c_{1}^{s,{avg}}} - \overset{\_}{c_{1}^{s,{avg}}}} \right\rbrack}} - {8D_{1}^{s}\overset{\_}{q_{1}^{avg}}}} = {- \overset{\_}{j_{1}}}$$\begin{matrix}{\frac{i_{app}}{{Fa}_{1}l_{1}} = {{k_{1}\left( \overset{\_}{c_{1}} \right)}_{1}^{\alpha_{c,1}}\left( {c_{1}^{s,\max} - \overset{\_}{c_{1}^{s,{surf}}}} \right)^{\alpha_{a,1}}\left( \overset{\_}{c_{1}^{s,{surf}}} \right)^{\alpha_{c,1}}\left( {{\exp\left( \frac{\alpha_{a,1}F\overset{\_}{\eta_{1}}}{RT} \right)} - {\exp\left( \frac{{- \alpha_{c,1}}F\overset{\_}{\eta_{1}}}{RT} \right)}} \right)}} \\{\overset{\_}{\eta_{1}} = {\overset{\_}{\phi_{s,1}} - \overset{\_}{\phi_{l,1}} - {U\left( \overset{\_}{c_{1}^{s,{surf}}} \right)}}}\end{matrix}$ Separator (Region 2)$\frac{d\overset{\_}{c_{2}}}{dt} = \frac{\frac{{- 2}{D\left( c_{12} \right)}\left( {\overset{\_}{c_{2}} - \overset{\_}{c_{1}}} \right)}{\frac{l_{1}}{\varepsilon_{1}^{b_{1}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} + \frac{2{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}}{\varepsilon_{2}l_{2}}$$\phi_{l,12} = {\left( \frac{{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}}\overset{\_}{\phi_{l,1}}} + {\frac{\varepsilon_{2}^{b_{2}}}{l_{2}}\overset{\_}{\phi_{l,2}}}}{\frac{\varepsilon_{1}^{b_{1}}}{l_{1}} + \frac{\varepsilon_{2}^{b_{2}}}{l_{2}}} \right) = 0}$Negative Electrode (Region 3)$\frac{d\overset{\_}{c_{3}}}{dt} = {\frac{\frac{{- 2}{D\left( c_{23} \right)}\left( {\overset{\_}{c_{3}} - \overset{\_}{c_{2}}} \right)}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}}}{\varepsilon_{3}l_{3}} - {\left( {1 - t_{+}^{0}} \right)\frac{i_{app}}{\varepsilon_{3}l_{3}}}}$$i_{app} = {{{- 2}{\kappa\left( c_{23} \right)}\left( \frac{\overset{\_}{\phi_{l,3}} - \overset{\_}{\phi_{l,2}}}{\frac{l_{3}}{\varepsilon_{3}^{b_{3}}} + \frac{l_{2}}{\varepsilon_{2}^{b_{2}}}} \right)} + {\frac{4R{T\left( {1 - t_{+}^{0}} \right)}}{F}{v\left( c_{23} \right)}{\kappa\left( c_{23} \right)}\frac{1}{c_{23}}\left( \frac{\overset{\_}{c_{3}} - \overset{\_}{c_{2}}}{\frac{l_{2}}{\varepsilon_{2}^{b_{2}}} + \frac{l_{3}}{\varepsilon_{3}^{b_{3}}}} \right)}}$$\frac{d\overset{\_}{c_{3}^{s,{avg}}}}{dt} = {{- 3}\frac{\overset{\_}{j_{3}}}{R_{3}}}$$\frac{d\overset{\_}{q_{3}^{avg}}}{dt} = {{{- 3}0\frac{D_{3}^{s}}{R_{3}^{2}}\overset{\_}{q_{3}^{avg}}} - {\frac{45}{2}\frac{\overset{\_}{j_{3}}}{R_{3}^{2}}}}$${{35{\frac{D_{3}^{s}}{R_{1}}\left\lbrack {\overset{\_}{c_{3}^{s,{surf}}} - \overset{\_}{c_{3}^{s,{avg}}}} \right\rbrack}} - {8D_{3}^{s}\overset{\_}{q_{3}^{avg}}}} = {- \overset{\_}{j_{3}}}$$\begin{matrix}{\frac{- i_{app}}{{Fa}_{3}l_{3}} = {{k_{3}\left( \overset{\_}{c_{3}} \right)}^{\alpha_{a,3}}\left( {c_{3}^{s,\max} - \overset{\_}{c_{3}^{s,{surf}}}} \right)^{\alpha_{a,3}}\left( \overset{\_}{c_{3}^{s,{surf}}} \right)^{\alpha_{c,3}}\left( {{\exp\left( \frac{\alpha_{a,3}F\overset{\_}{\eta_{3}}}{RT} \right)} - {\exp\left( \frac{{- \alpha_{c,3}}F\overset{\_}{\eta_{3}}}{RT} \right)}} \right)}} \\{\overset{\_}{\eta_{3}} = {\overset{\_}{\phi_{s,3}} - \overset{\_}{\phi_{l,3}} - {U\left( \overset{\_}{c_{3}^{s,{surf}}} \right)}}}\end{matrix}$

TABLE V Representative computational performance metrics for the TankModel. Model and Number of Computation Implementation DAEs time (ms)*p2D-Finite Difference 986  1493 (50, 35, 50) p2D-Reformulated 30-1505-20 Tank Model  15** 2.1 *Average of N = 3 simulation runs for a 1 Cdischarge. It must be noted that the simulation time and memoryconsumption is a strong function of the computing environment as well aserror tolerances for both initialization and simulation. **Including theequation for cell voltage

List of Symbols Dependent Variables C Electrolyte Concentration c^(s)Solid Phase Concentration ϕ_(l) Liquid Phase Potential ϕ_(s) Solid PhasePotential V_(cell) Cell Voltage j Pore-wall flux N Electrolyte molarflux Other Superscripts surf Pertaining to the surface of the particlein the solid phase — Pertaining to the average over the volume of aporous domain s, avg Pertaining to the average over the volume of thesolid particle Other Subscripts i Pertaining to region i where i ϵ {1,2, 3} ij Pertaining to the interface between regions i and j, where i, jϵ {1, 2, 3} i, ij Pertaining to the interface between regions i and j onthe side of region i where i, j ϵ {1, 2, 3}

The invention claimed is:
 1. A battery management system, comprising: aconnector for electrically coupling a battery to the battery managementsystem; at least one sensor configured to detect a battery state; aprogrammable chip configured to control at least one of charging anddischarging of the battery; and a controller device configured to:receive at least one battery state from the at least one sensor; providethe at least one battery state as input to a tanks-in-series model thatrepresents the battery; and provide at least one output of thetanks-in-series model to the programmable chip for controlling at leastone of charging and discharging of the battery.
 2. The batterymanagement system of claim 1, wherein the tanks-in-series modelincludes: a tank model that represents an anode of the battery; a tankmodel that represents a separator of the battery; and a tank model thatrepresents a cathode of the battery.
 3. The battery management system ofclaim 2, wherein each tank model includes a volume-averagedrepresentation of a P2D model of a corresponding domain.
 4. The batterymanagement system of claim 3, wherein the tanks-in-series model isderived by: integrating a pseudo 2-dimensional (P2D) modelrepresentation over a volume of each domain; and approximatinginterfacial fluxes between the domains, wherein a corresponding gradientis expressed as an average value minus an interfacial value of adependent variable over a length scale.
 5. The battery management systemof claim 4, wherein the derivation of the tanks-in-series model furtherincludes determining the interfacial values in a manner that preservesmass and charge conservation between each tank model.
 6. The batterymanagement system of claim 4, wherein integrating the P2D modelrepresentation over the volume of each domain results in a set ofvolume-averaged quantities with temporal evolution of each averagedvariable expressed as an overall balance comprising approximatedinternal source terms and interfacial fluxes.
 7. The battery managementsystem of claim 1, wherein the at least one battery state includes atleast one of a terminal voltage of the battery and a surface temperatureof the battery.
 8. The battery management system of claim 1, wherein theoutput includes at least one of a state of charge and a state of healthof the battery.
 9. A method of managing at least one of charging anddischarging of a battery, the method comprising: receiving, by acontroller device, at least one battery state of the battery from atleast one sensor; providing, by the controller device, the at least onebattery state as input to a tanks-in-series model that represents thebattery; and controlling, by the controller device, at least one ofcharging or discharging of the battery by using at least one output ofthe tanks-in-series model to determine at least one of a current or avoltage for use in the at least one of charging or discharging of thebattery; wherein the tanks-in-series model includes: a tank model thatrepresents an anode of the battery; a tank model that represents aseparator of the battery; and a tank model that represents a cathode ofthe battery.
 10. The method of claim 9, wherein the tank model thatrepresents the anode of the battery includes a volume-averagedrepresentation of a concentration of lithium ions and potential in anelectrolyte of the anode; the tank model that represents the separatorof the battery includes a volume-averaged representation of aconcentration of lithium ions and potential in an electrolyte of theseparator; and the tank model that represents the cathode of the batteryincludes a volume-averaged representation of a concentration of lithiumions and potential in an electrolyte of the cathode.
 11. The method ofclaim 10, wherein each tank model includes a volume-averagedrepresentation of a pseudo 2-dimensional (p2D) model of a correspondingdomain.
 12. The method of claim 11, wherein the tanks-in-series model isderived by: integrating a p2D model representation over a volume of eachdomain; and approximating interfacial fluxes between the domains,wherein a corresponding gradient is expressed as an average value minusan interfacial value of a dependent variable over a length scale. 13.The method of claim 12, wherein the derivation of the tanks-in-seriesmodel further includes determining the interfacial values in a mannerthat preserves mass and charge conservation between each tank model. 14.The method of claim 12, wherein integrating the p2D model representationover the volume of each domain results in a set of volume-averagedquantities with temporal evolution of each averaged variable expressedas an overall balance comprising approximated internal source terms andinterfacial fluxes.
 15. The method of claim 9, wherein the at least onebattery state includes at least one of a terminal voltage of the batteryand a surface temperature of the battery.
 16. The method of claim 9,wherein the output includes at least one of a state of charge and astate of health of the battery.
 17. A computer-implemented method ofderiving a model for controlling at least one of charging anddischarging a battery, the method comprising: integrating, by acomputing device, pseudo 2-dimensional (p2D) model representations overvolumes of a cathode domain of the battery, an anode domain of thebattery, and a separator domain of the battery; approximating, by thecomputing device, interfacial fluxes between the domains, wherein acorresponding gradient is expressed as an average value minus aninterfacial value of a dependent variable over a length scale to createa tanks-in-series representation of the battery; and providing, by thecomputing device, the tanks-in-series representation of the battery asthe model for controlling at least one of charging and discharging ofthe battery.
 18. The computer-implemented method of claim 17, whereincreating the tanks-in-series representation of the battery includesdetermining, by the computing device, the interfacial values in a mannerthat preserves mass and charge conservation between each domain of thebattery.
 19. The computer-implemented method of claim 17, whereinintegrating the p2D model representations over the volumes of thedomains results in a set of volume-averaged quantities with temporalevolution of each averaged variable expressed as an overall balancecomprising approximated internal source terms and interfacial fluxes.20. The computer-implemented method of claim 17, wherein the modelaccepts at least one of a terminal voltage of the battery and a surfacetemperature of the battery as input, and wherein the model outputs atleast one of a state of charge of the battery and a health of thebattery.